*$ CREATE MGDRAW.FOR
*COPY MGDRAW
*                                                                      *
*=== mgdraw ===========================================================*
*                                                                      *
      SUBROUTINE MGDRAW ( ICODE, MREG )

      INCLUDE '(DBLPRC)'
      INCLUDE '(DIMPAR)'
      INCLUDE '(IOUNIT)'
*
*----------------------------------------------------------------------*
*                                                                      *
*     Copyright (C) 1990-2006      by        Alfredo Ferrari           *
*     All Rights Reserved.                                             *
*                                                                      *
*                                                                      *
*     MaGnetic field trajectory DRAWing: actually this entry manages   *
*                                        all trajectory dumping for    *
*                                        drawing                       *
*                                                                      *
*     Created on   01 march 1990   by        Alfredo Ferrari           *
*                                              INFN - Milan            *
*     Last change  05-may-06       by        Alfredo Ferrari           *
*                                              INFN - Milan            *
*                                                                      *
*----------------------------------------------------------------------*
*
      INCLUDE '(CASLIM)'
      INCLUDE '(COMPUT)'
      INCLUDE '(SOURCM)'
      INCLUDE '(FHEAVY)'
      INCLUDE '(FLKSTK)'
      INCLUDE '(GENSTK)'
      INCLUDE '(MGDDCM)'
      INCLUDE '(PAPROP)'
      INCLUDE '(QUEMGD)'
      INCLUDE '(SUMCOU)'
      INCLUDE '(TRACKR)'
      INCLUDE '(LTCLCM)'
*
      DIMENSION DTQUEN ( MXTRCK, MAXQMG )
*
      CHARACTER*20 FILNAM
      INTEGER I
      LOGICAL LFCOPE
      SAVE LFCOPE
      DATA LFCOPE / .FALSE. /
      CHARACTER*8 NAMEIN
      CHARACTER*8 NAMEOUT
      CHARACTER*8 NAMELAT
*
*----------------------------------------------------------------------*
*                                                                      *
*     Icode = 1: call from Kaskad                                      *
*     Icode = 2: call from Emfsco                                      *
*     Icode = 3: call from Kasneu                                      *
*     Icode = 4: call from Kashea                                      *
*     Icode = 5: call from Kasoph                                      *
*                                                                      *
*----------------------------------------------------------------------*
*                                                                      *
      IF ( .NOT. LFCOPE ) THEN
         LFCOPE = .TRUE.
         IF ( KOMPUT .EQ. 2 ) THEN
            FILNAM = '/'//CFDRAW(1:8)//' DUMP A'
         ELSE
            FILNAM = CFDRAW
         END IF
         OPEN ( UNIT = IODRAW, FILE = FILNAM, STATUS = 'NEW', FORM =
     &          'UNFORMATTED' )
      END IF
*      WRITE (IODRAW) NTRACK, MTRACK, JTRACK, SNGL (ETRACK),
*     &               SNGL (WTRACK)
*      WRITE (IODRAW) ( SNGL (XTRACK (I)), SNGL (YTRACK (I)),
*     &                 SNGL (ZTRACK (I)), I = 0, NTRACK ),
*     &               ( SNGL (DTRACK (I)), I = 1, MTRACK ),
*     &                 SNGL (CTRACK)
*       PRINT *,'NTrSeg=',NTRACK,' PartId=',JTRACK
*       DO I=0,NTRACK
*       PRINT *,'X=',XTRACK(I),' Y=',YTRACK(I),' Z=',ZTRACK(I)
*       END DO
*  +-------------------------------------------------------------------*
*  |  Quenching is activated
      IF ( LQEMGD ) THEN
         IF ( MTRACK .GT. 0 ) THEN
            RULLL  = ZERZER
            CALL QUENMG ( ICODE, MREG, RULLL, DTQUEN )
*           WRITE (IODRAW) ( ( SNGL (DTQUEN (I,JBK)), I = 1, MTRACK ),
*    &                         JBK = 1, NQEMGD )
         END IF
      END IF
*  |  End of quenching
*  +-------------------------------------------------------------------*
      RETURN
*
*======================================================================*
*                                                                      *
*     Boundary-(X)crossing DRAWing:                                    *
*                                                                      *
*     Icode = 1x: call from Kaskad                                     *
*             19: boundary crossing                                    *
*     Icode = 2x: call from Emfsco                                     *
*             29: boundary crossing                                    *
*     Icode = 3x: call from Kasneu                                     *
*             39: boundary crossing                                    *
*     Icode = 4x: call from Kashea                                     *
*             49: boundary crossing                                    *
*     Icode = 5x: call from Kasoph                                     *
*             59: boundary crossing                                    *
*                                                                      *
*======================================================================*
*                                                                      *
      ENTRY BXDRAW ( ICODE, MREG, NEWREG, XSCO, YSCO, ZSCO )
      IF (ISPUSR(MKBMX2).EQ.1) THEN
        CALL GEOR2N(MREG,NAMEIN,ERR)
        CALL GEOR2N(NEWREG,NAMEOUT,ERR)
        CALL GEOL2N(MLATTC,NAMELAT,IRTLAT,ERR)
        PRINT *,"At event ",Ncase," particle moving form region ",NAMEIN
     +       ," to region ",NAMEOUT
        IF ((MOD(MREG,4).EQ.1).AND.(MOD(NEWREG,4).EQ.3)) THEN
            PRINT *,"Magnet entering  ",NAMELAT,
     +      "at (x,y,z):",XSCO,YSCO,ZSCO
        ELSE IF ((MOD(MREG,4).EQ.3).AND.(MOD(NEWREG,4).EQ.1)) THEN
            PRINT *,"Magnet exiting  ",NAMELAT,
     +      "at (x,y,z):",XSCO,YSCO,ZSCO
        END IF	
      END IF 
*      IF ((JTRACK.EQ.1).AND.(ISPUSR(MKBMX2).EQ.1)) THEN
*      PRINT *,MREG,"-->",NEWREG," at LATTICENo=",MLATTC
*      PRINT *,"EventNo=",Ncase," at position(x,y,x)",XSCO,YSCO,ZSCO
*      END IF
*     IF ((JTRACK.EQ.1).AND.(ISPUSR(MKBMX2).EQ.1)) THEN
*	IF ((MREG.EQ.1).AND.(NEWREG.EQ.3).AND.(MLATTC.EQ.1)) THEN
*	PRINT *,"B3Ta input at(x,y,z)=",XSCO,YSCO,ZSCO
*	END IF
*	IF ((MREG.EQ.3).AND.(NEWREG.EQ.1).AND.(MLATTC.EQ.1)) THEN
*	PRINT *,"B3Ta output at(x,y,z)=",XSCO,YSCO,ZSCO
*	END IF
*	IF ((MREG.EQ.1).AND.(NEWREG.EQ.3).AND.(MLATTC.EQ.2)) THEN
*	PRINT *,"B3Tb input at(x,y,z)=",XSCO,YSCO,ZSCO
*	END IF
*	IF ((MREG.EQ.3).AND.(NEWREG.EQ.1).AND.(MLATTC.EQ.2)) THEN
*	PRINT *,"B3Tb output at(x,y,z)=",XSCO,YSCO,ZSCO
*	END IF
*	IF ((MREG.EQ.5).AND.(NEWREG.EQ.7).AND.(MLATTC.EQ.3)) THEN
*	PRINT *,"B1 input at(x,y,z)=",XSCO,YSCO,ZSCO
*	END IF
*	IF ((MREG.EQ.7).AND.(NEWREG.EQ.5).AND.(MLATTC.EQ.3)) THEN
*	PRINT *,"B1 output at(x,y,z)=",XSCO,YSCO,ZSCO
*	END IF
*	IF ((MREG.EQ.5).AND.(NEWREG.EQ.7).AND.(MLATTC.EQ.4)) THEN
*	PRINT *,"B11 input at(x,y,z)=",XSCO,YSCO,ZSCO
*	END IF
*	IF ((MREG.EQ.7).AND.(NEWREG.EQ.5).AND.(MLATTC.EQ.4)) THEN
*	PRINT *,"B11 output at(x,y,z)=",XSCO,YSCO,ZSCO
*	END IF
*	IF ((MREG.EQ.17).AND.(NEWREG.EQ.19).AND.(MLATTC.EQ.5)) THEN
*	PRINT *,"Q2 input at(x,y,z)=",XSCO,YSCO,ZSCO
*	END IF
*	IF ((MREG.EQ.19).AND.(NEWREG.EQ.17).AND.(MLATTC.EQ.5)) THEN
*	PRINT *,"Q2 output at(x,y,z)=",XSCO,YSCO,ZSCO
*	END IF
*	IF ((MREG.EQ.25).AND.(NEWREG.EQ.27).AND.(MLATTC.EQ.6)) THEN
*	PRINT *,"Q1 input at(x,y,z)=",XSCO,YSCO,ZSCO
*	END IF
*	IF ((MREG.EQ.27).AND.(NEWREG.EQ.25).AND.(MLATTC.EQ.6)) THEN
*	PRINT *,"Q1 output at(x,y,z)=",XSCO,YSCO,ZSCO
*	END IF
*	IF ((MREG.EQ.29).AND.(NEWREG.EQ.31).AND.(MLATTC.EQ.7)) THEN
*	PRINT *,"C1 input at(x,y,z)=",XSCO,YSCO,ZSCO
*	END IF
*	IF ((MREG.EQ.31).AND.(NEWREG.EQ.29).AND.(MLATTC.EQ.7)) THEN
*	PRINT *,"C1 output at(x,y,z)=",XSCO,YSCO,ZSCO
*	END IF
*	IF ((MREG.EQ.37).AND.(NEWREG.EQ.39).AND.(MLATTC.EQ.8)) THEN
*	PRINT *,"C2 input at(x,y,z)=",XSCO,YSCO,ZSCO
*	END IF
*	IF ((MREG.EQ.39).AND.(NEWREG.EQ.37).AND.(MLATTC.EQ.8)) THEN
*	PRINT *,"C2 output at(x,y,z)=",XSCO,YSCO,ZSCO
*	END IF
*	IF ((MREG.EQ.17).AND.(NEWREG.EQ.19).AND.(MLATTC.EQ.9)) THEN
*	PRINT *,"Q3 input at(x,y,z)=",XSCO,YSCO,ZSCO
*	END IF
*	IF ((MREG.EQ.19).AND.(NEWREG.EQ.17).AND.(MLATTC.EQ.9)) THEN
*	PRINT *,"Q3 output at(x,y,z)=",XSCO,YSCO,ZSCO
*	END IF
*	IF ((MREG.EQ.13).AND.(NEWREG.EQ.15).AND.(MLATTC.EQ.10)) THEN
*	PRINT *,"B2a input at(x,y,z)=",XSCO,YSCO,ZSCO
*	END IF
*	IF ((MREG.EQ.15).AND.(NEWREG.EQ.13).AND.(MLATTC.EQ.10)) THEN
*	PRINT *,"B2a output at(x,y,z)=",XSCO,YSCO,ZSCO
*	END IF
*	IF ((MREG.EQ.13).AND.(NEWREG.EQ.15).AND.(MLATTC.EQ.11)) THEN
*	PRINT *,"B2b input at(x,y,z)=",XSCO,YSCO,ZSCO
*	END IF
*	IF ((MREG.EQ.15).AND.(NEWREG.EQ.13).AND.(MLATTC.EQ.11)) THEN
*	PRINT *,"B2b output at(x,y,z)=",XSCO,YSCO,ZSCO
*	END IF
*	IF ((MREG.EQ.13).AND.(NEWREG.EQ.15).AND.(MLATTC.EQ.12)) THEN
*	PRINT *,"B2c input at(x,y,z)=",XSCO,YSCO,ZSCO
*	END IF
*	IF ((MREG.EQ.15).AND.(NEWREG.EQ.13).AND.(MLATTC.EQ.12)) THEN
*	PRINT *,"B2c output at(x,y,z)=",XSCO,YSCO,ZSCO
*	END IF
*	IF ((MREG.EQ.49).AND.(NEWREG.EQ.51).AND.(MLATTC.EQ.13)) THEN
*	PRINT *,"C3.C4 input at(x,y,z)=",XSCO,YSCO,ZSCO
*	END IF
*	IF ((MREG.EQ.51).AND.(NEWREG.EQ.49).AND.(MLATTC.EQ.13)) THEN
*	PRINT *,"C3,C4 output at(x,y,z)=",XSCO,YSCO,ZSCO
*	END IF
*	IF ((MREG.EQ.13).AND.(NEWREG.EQ.15).AND.(MLATTC.EQ.14)) THEN
*	PRINT *,"B3a input at(x,y,z)=",XSCO,YSCO,ZSCO
*	END IF
*	IF ((MREG.EQ.15).AND.(NEWREG.EQ.13).AND.(MLATTC.EQ.14)) THEN
*	PRINT *,"B3a output at(x,y,z)=",XSCO,YSCO,ZSCO
*	END IF
*	IF ((MREG.EQ.13).AND.(NEWREG.EQ.15).AND.(MLATTC.EQ.15)) THEN
*	PRINT *,"B3b input at(x,y,z)=",XSCO,YSCO,ZSCO
*	END IF
*	IF ((MREG.EQ.15).AND.(NEWREG.EQ.13).AND.(MLATTC.EQ.15)) THEN
*	PRINT *,"B3b output at(x,y,z)=",XSCO,YSCO,ZSCO
*	END IF
*	IF ((MREG.EQ.13).AND.(NEWREG.EQ.15).AND.(MLATTC.EQ.16)) THEN
*	PRINT *,"B3c input at(x,y,z)=",XSCO,YSCO,ZSCO
*	END IF
*	IF ((MREG.EQ.15).AND.(NEWREG.EQ.13).AND.(MLATTC.EQ.16)) THEN
*	PRINT *,"B3c output at(x,y,z)=",XSCO,YSCO,ZSCO
*	END IF
*	IF ((MREG.EQ.17).AND.(NEWREG.EQ.19).AND.(MLATTC.EQ.17)) THEN
*	PRINT *,"Q4 input at(x,y,z)=",XSCO,YSCO,ZSCO
*	END IF
*	IF ((MREG.EQ.19).AND.(NEWREG.EQ.17).AND.(MLATTC.EQ.17)) THEN
*	PRINT *,"Q4 output at(x,y,z)=",XSCO,YSCO,ZSCO
*	END IF
*	IF ((MREG.EQ.17).AND.(NEWREG.EQ.19).AND.(MLATTC.EQ.18)) THEN
*	PRINT *,"Q5 input at(x,y,z)=",XSCO,YSCO,ZSCO
*	END IF
*	IF ((MREG.EQ.19).AND.(NEWREG.EQ.17).AND.(MLATTC.EQ.18)) THEN
*	PRINT *,"Q5 output at(x,y,z)=",XSCO,YSCO,ZSCO
*	END IF
*	IF ((MREG.EQ.53).AND.(NEWREG.EQ.55).AND.(MLATTC.EQ.19)) THEN
*	PRINT *,"C5,C6 input at(x,y,z)=",XSCO,YSCO,ZSCO
*	END IF
*	IF ((MREG.EQ.55).AND.(NEWREG.EQ.53).AND.(MLATTC.EQ.19)) THEN
*	PRINT *,"C5,C6 output at(x,y,z)=",XSCO,YSCO,ZSCO
*	END IF
*	IF ((MREG.EQ.21).AND.(NEWREG.EQ.23).AND.(MLATTC.EQ.20)) THEN
*	PRINT *,"Q6a input at(x,y,z)=",XSCO,YSCO,ZSCO
*	END IF
*	IF ((MREG.EQ.23).AND.(NEWREG.EQ.21).AND.(MLATTC.EQ.20)) THEN
*	PRINT *,"Q6a output at(x,y,z)=",XSCO,YSCO,ZSCO
*	END IF
*     END IF
*      PRINT *,"Ooppss Event=",Ncase," LAT=",MLATTC,"REG=",MREG,NEWREG
      RETURN
*
*======================================================================*
*                                                                      *
*     Event End DRAWing:                                               *
*                                                                      *
*======================================================================*
*                                                                      *
      ENTRY EEDRAW ( ICODE )
      RETURN
*
*======================================================================*
*                                                                      *
*     ENergy deposition DRAWing:                                       *
*                                                                      *
*     Icode = 1x: call from Kaskad                                     *
*             10: elastic interaction recoil                           *
*             11: inelastic interaction recoil                         *
*             12: stopping particle                                    *
*             13: pseudo-neutron deposition                            *
*             14: escape                                               *
*             15: time kill                                            *
*     Icode = 2x: call from Emfsco                                     *
*             20: local energy deposition (i.e. photoelectric)         *
*             21: below threshold, iarg=1                              *
*             22: below threshold, iarg=2                              *
*             23: escape                                               *
*             24: time kill                                            *
*     Icode = 3x: call from Kasneu                                     *
*             30: target recoil                                        *
*             31: below threshold                                      *
*             32: escape                                               *
*             33: time kill                                            *
*     Icode = 4x: call from Kashea                                     *
*             40: escape                                               *
*             41: time kill                                            *
*             42: delta ray stack overflow                             *
*     Icode = 5x: call from Kasoph                                     *
*             50: optical photon absorption                            *
*             51: escape                                               *
*             52: time kill                                            *
*                                                                      *
*======================================================================*
*                                                                      *
      ENTRY ENDRAW ( ICODE, MREG, RULL, XSCO, YSCO, ZSCO )
      IF ( .NOT. LFCOPE ) THEN
         LFCOPE = .TRUE.
         IF ( KOMPUT .EQ. 2 ) THEN
            FILNAM = '/'//CFDRAW(1:8)//' DUMP A'
         ELSE
            FILNAM = CFDRAW
         END IF
         OPEN ( UNIT = IODRAW, FILE = FILNAM, STATUS = 'NEW', FORM =
     &          'UNFORMATTED' )
      END IF
*      WRITE (IODRAW)  0, ICODE, JTRACK, SNGL (ETRACK), SNGL (WTRACK)
*      WRITE (IODRAW)  SNGL (XSCO), SNGL (YSCO), SNGL (ZSCO), SNGL (RULL)
*  +-------------------------------------------------------------------*
*  |  Quenching is activated : calculate quenching factor
*  |  and store quenched energy in DTQUEN(1, jbk)
      IF ( LQEMGD ) THEN
         RULLL = RULL
         CALL QUENMG ( ICODE, MREG, RULLL, DTQUEN )
*         WRITE (IODRAW) ( SNGL (DTQUEN(1, JBK)), JBK = 1, NQEMGD )
      END IF
*  |  end quenching
*  +-------------------------------------------------------------------*
      RETURN
*
*======================================================================*
*                                                                      *
*     SOurce particle DRAWing:                                         *
*                                                                      *
*======================================================================*
*
      ENTRY SODRAW
      IF ( .NOT. LFCOPE ) THEN
         LFCOPE = .TRUE.
         IF ( KOMPUT .EQ. 2 ) THEN
            FILNAM = '/'//CFDRAW(1:8)//' DUMP A'
         ELSE
            FILNAM = CFDRAW
         END IF
         OPEN ( UNIT = IODRAW, FILE = FILNAM, STATUS = 'NEW', FORM =
     &          'UNFORMATTED' )
      END IF
*     WRITE (IODRAW) -NCASE, NPFLKA, NSTMAX, SNGL (TKESUM),
*    &                SNGL (WEIPRI)
*  +-------------------------------------------------------------------*
*  |  (Radioactive) isotope: it works only for 1 source particle on
*  |  the stack for the time being
      IF ( ILOFLK (NPFLKA) .GE. 100000 .AND. LRADDC (NPFLKA) ) THEN
         IARES  = MOD ( ILOFLK (NPFLKA), 100000  )  / 100
         IZRES  = MOD ( ILOFLK (NPFLKA), 10000000 ) / 100000
         IISRES = ILOFLK (NPFLKA) / 10000000
         IONID  = ILOFLK (NPFLKA)
*         WRITE (IODRAW) ( IONID,SNGL(-TKEFLK(I)),
*     &                    SNGL (WTFLK(I)), SNGL (XFLK (I)),
*     &                    SNGL (YFLK (I)), SNGL (ZFLK (I)),
*     &                    SNGL (TXFLK(I)), SNGL (TYFLK(I)),
*     &                    SNGL (TZFLK(I)), I = 1, NPFLKA )
*  |
*  +-------------------------------------------------------------------*
*  |  Patch for heavy ions: it works only for 1 source particle on
*  |  the stack for the time being
      ELSE IF ( ABS (ILOFLK (NPFLKA)) .GE. 10000 ) THEN
         IONID = ILOFLK (NPFLKA)
         CALL DCDION ( IONID )
*         WRITE (IODRAW) ( IONID,SNGL(TKEFLK(I)+AMNHEA(-IONID)),
*     &                    SNGL (WTFLK(I)), SNGL (XFLK (I)),
*     &                    SNGL (YFLK (I)), SNGL (ZFLK (I)),
*     &                    SNGL (TXFLK(I)), SNGL (TYFLK(I)),
*     &                    SNGL (TZFLK(I)), I = 1, NPFLKA )
*  |
*  +-------------------------------------------------------------------*
*  |  Patch for heavy ions: ???
*     ELSE IF ( ILOFLK (NPFLKA) .LT. -6 ) THEN
*         WRITE (IODRAW) ( IONID,SNGL(TKEFLK(I)+AMNHEA(-ILOFLK(NPFLKA))),
*     &                    SNGL (WTFLK(I)), SNGL (XFLK (I)),
*     &                    SNGL (YFLK (I)), SNGL (ZFLK (I)),
*     &                    SNGL (TXFLK(I)), SNGL (TYFLK(I)),
*     &                    SNGL (TZFLK(I)), I = 1, NPFLKA )
*  |
*  +-------------------------------------------------------------------*
*  |
*      ELSE
*         WRITE (IODRAW) ( ILOFLK(I), SNGL (TKEFLK(I)+AM(ILOFLK(I))),
*     &                    SNGL (WTFLK(I)), SNGL (XFLK (I)),
*     &                    SNGL (YFLK (I)), SNGL (ZFLK (I)),
*     &                    SNGL (TXFLK(I)), SNGL (TYFLK(I)),
*     &                    SNGL (TZFLK(I)), I = 1, NPFLKA )
      END IF
*  |
*  +-------------------------------------------------------------------*
      RETURN
*
*======================================================================*
*                                                                      *
*     USer dependent DRAWing:                                          *
*                                                                      *
*     Icode = 10x: call from Kaskad                                    *
*             100: elastic   interaction secondaries                   *
*             101: inelastic interaction secondaries                   *
*             102: particle decay  secondaries                         *
*             103: delta ray  generation secondaries                   *
*             104: pair production secondaries                         *
*             105: bremsstrahlung  secondaries                         *
*             110: decay products                                      *
*     Icode = 20x: call from Emfsco                                    *
*             208: bremsstrahlung secondaries                          *
*             210: Moller secondaries                                  *
*             212: Bhabha secondaries                                  *
*             214: in-flight annihilation secondaries                  *
*             215: annihilation at rest   secondaries                  *
*             217: pair production        secondaries                  *
*             219: Compton scattering     secondaries                  *
*             221: photoelectric          secondaries                  *
*             225: Rayleigh scattering    secondaries                  *
*     Icode = 30x: call from Kasneu                                    *
*             300: interaction secondaries                             *
*     Icode = 40x: call from Kashea                                    *
*             400: delta ray  generation secondaries                   *
*  For all interactions secondaries are put on GENSTK common (kp=1,np) *
*  but for KASHEA delta ray generation where only the secondary elec-  *
*  tron is present and stacked on FLKSTK common for kp=npflka          *
*                                                                      *
*======================================================================*
*
      ENTRY USDRAW ( ICODE, MREG, XSCO, YSCO, ZSCO )
      IF ( .NOT. LFCOPE ) THEN
         LFCOPE = .TRUE.
         IF ( KOMPUT .EQ. 2 ) THEN
            FILNAM = '/'//CFDRAW(1:8)//' DUMP A'
         ELSE
            FILNAM = CFDRAW
         END IF
         OPEN ( UNIT = IODRAW, FILE = FILNAM, STATUS = 'NEW', FORM =
     &          'UNFORMATTED' )
      END IF
* No output by default:
      RETURN
*=== End of subrutine Mgdraw ==========================================*
      END

