*
* $Id$
*
* $Log$
*
*CMZ :  2.01/00 06/03/2000  13.07.03  by  Francesco Cafagna
*CMZ :  2.00/00 03/03/2000  15.39.06  by  Francesco Cafagna
*CMZ :  1.02/00 22/10/99  11.10.56  by  Francesco Cafagna
*CMZ :  1.01/00 21/05/96  08.53.57  by  Francesco Cafagna
*-- Author :    Paolo (The Magic) Papini   23/04/96
#if defined(GPAMELA_PAOLO)
      SUBROUTINE GPSPEHD(NHIT,NPLANE,XI,YI,ZI,XO,YO,ZO,PATH,DE,LKDN)
#endif
#if !defined(GPAMELA_PAOLO)
      SUBROUTINE GPSPEHD(NHIT,NUMV,NLEV,XI,YI,ZI,XO,YO,ZO,PATH,DE,
     +                   LKDN)
#endif
************************************************************************
*                                                                      *
* To perform hit/digit conversion for spectrometer.                    *
* This routine fills the LKD1,LKD2,LKD3,LKD4 array                     *
* and it up-date LKND as total number of digit                         *
*                                                                      *
* Variables definition:                                                *
* IN:                                                                  *
*     NHIT - hit number                                                *
*     NPLANE - plane number                                            *
*     XI,YI,ZI,XO,YO,ZO - input and output hit coordinates             *
*     PATH,DE - hit path lenght and dE/dx                              *
*                                                                      *
* OUT:                                                                 *
*     LKDN - total number of digit                                     *
*                                                                      *
* Called by: GPDSPE                                                    *
* Author: Paolo (The Magic) Papini, 23/04/96 15.45.26                  *
*                                                                      *
************************************************************************
#include "gcunit.inc"
#include "gpphit.inc"
#include "gppdigi.inc"
      PARAMETER (NH=NHSPE,NHM=NHMSPE,NVOL=4)
#include "gpspe.inc"
#include "gpsed.inc"
#include "gpgeo.inc"
#include "gpdspe.inc"
*
* Data block for digitalization
*
      INTEGER NHIT,NPLANE,LKDN,NSIL,NSTRIP,NADC,NEWD
      REAL XSIL,YSIL,XBORD,YBORD
      REAL XI,YI,XO,YO,PATH,DE,CH_DIG(2048)
      REAL RR(3),BB(3)
#if !defined(GPAMELA_PAOLO)
      INTEGER LNAM,LNUM
      DIMENSION NUMV(NLEV),LNAM(NVOL),LNUM(NVOL)
      REAL XYZISI(3),XYZOSI(3),XYZIMA(3),XYZOMA(3),XYZMO(3),XYZSO(3)
#endif
*
      BB(1)=0.
      BB(2)=0.
      BB(3)=0.
#if defined(GPAMELA_MAGN)
      RR(1)=(XI+XO)/2.
      RR(2)=(YI+YO)/2.
      RR(3)=(ZI+ZO)/2.
      CALL GUFLD(RR,BB)
#endif
      NSEPX=1+INT(ABS(XO-XI)/STEPX)
      NSEPY=1+INT(ABS(YO-YI)/STEPY)
      NSEPZ=1+INT(ABS(ZO-ZI)/STEPZ)
#if defined(GPAMELA_PAOLO)
C*
C* X-view digitalization
C*
C      NSIL=INT(((XI+XO)/2.+3.*TSPA(1))/(2.*TSPA(1)))+1
C* go way if no in the silicon
C      IF(NSIL.LT.1.OR.NSIL.GT.3) GOTO 10
C* n.division of hit path
C      NSEP=MAX0(NSEPX,NSEPZ)
C* energy loss in each path division
C      DDE=DE/FLOAT(NSEP)
C* XI-ccor. in the silicon
C      XSIL_I=MOD(XI+3.*TSPA(1),2.*TSPA(1))
C* XO-coor. in the silicon
C      XSIL_O=MOD(XO+3.*TSPA(1),2.*TSPA(1))
C* x-coor. of sensitive edge
C      XBORD=(2.*TSPA(1)-NXSTR*WXSTR)/2.
C* XI of the strip
C      NSTR_I=INT((XSIL_I-XBORD)/WXSTR)+1
C* XO of the strip
C      NSTR_O=INT((XSIL_O-XBORD)/WXSTR)+1
C      NSTR_MIN=MIN0(NSTR_I,NSTR_O)-10
C* minimum strip to be calculated
C      IF(NSTR_MIN.LT.1) NSTR_MIN=1
C      NSTR_MAX=MAX0(NSTR_I,NSTR_O)+10
C* maximum strip to be calculated
C      IF(NSTR_MAX.GT.NXSTR) NSTR_MAX=NXSTR
C+SELF
#endif
#if !defined(GPAMELA_PAOLO)
*
* Fill up the GCVOLU common block
*
*      CALL GFPATH (ISSPE,IDTSPA,NUMV,NLEV,LNAM,LNUM)
c      LNAM(1) =4HPAME
c      LNAM(2) =4HSPEB
      LNAM(1) =4HTRPB
      LNAM(2) =4HTRPL
      LNAM(3) =4HTRSL
      LNAM(4) =4HTSPA
      CALL GLVOLU (NLEV,LNAM,NUMV,IER)
*
* X-view digitalization
*
      NSIL=NUMV(NLEV-1)
* go way if no in the silicon
      IF(NSIL.LT.1.OR.NSIL.GT.3) GOTO 10
*
* Calculate the stripes
*
      NPLANE = NUMV(NLEV-3)
* n.division of hit path
      NSEP=MAX0(NSEPX,NSEPZ)
* energy loss in each path division
      DDE=DE/FLOAT(NSEP)
* trasform coordinates from MARS to silicon center
      XYZIMA(1) = XI
      XYZIMA(2) = YI
      XYZIMA(3) = ZI
      CALL GMTOD(XYZIMA,XYZISI,1)
      XYZOMA(1) = XO
      XYZOMA(2) = YO
      XYZOMA(3) = ZO
      CALL GMTOD(XYZOMA,XYZOSI,1)
* xy-coor. of sensitive edge
      XBORD=(2.*TSPA(1)-NXSTR*WXSTR)/2.
      YBORD=(2.*TSPA(2)-NYSTR*WYSTR)/2.
* xyz-coord. of the paddle origin in MARS
      CALL VZERO(XYZSO,3)
      CALL GDTOM(XYZSO,XYZMO,1)
* Move origin of coordinates to the paddle edge, a guard zone added
      XYZISI(1) = XYZISI(1) + (TSPA(1) - XBORD)
      XYZISI(2) = XYZISI(2) + (TSPA(2) - YBORD)
*
      XYZOSI(1) = XYZOSI(1) + (TSPA(1) - XBORD)
      XYZOSI(2) = XYZOSI(2) + (TSPA(2) - YBORD)
* Calculated the strip number now
      NSTR_I=INT(XYZISI(1)/WXSTR)+1
      NSTR_O=INT(XYZOSI(1)/WXSTR)+1
      NSTR_MIN=MIN0(NSTR_I,NSTR_O)-10
* minimum strip to be calculated
      IF(NSTR_MIN.LT.1) NSTR_MIN=1
      NSTR_MAX=MAX0(NSTR_I,NSTR_O)+10
* maximum strip to be calculated
      IF(NSTR_MAX.GT.NXSTR) NSTR_MAX=NXSTR
#endif
C      PRINT *,'X DIGIT, N.HIT:',NHIT
C      PRINT *,'NSTR_MIN,NSTR_MAX,NSEP',NSTR_MIN,NSTR_MAX,NSEP

* loop on the strip (for the diffusion of Q)
      DO IS=NSTR_MIN,NSTR_MAX
* x-center of the strip IS
         XSTR=XBORD+(IS-0.5)*WXSTR+(2*NSIL-5)*TSPA(1)
C+SELF,IF=PAOLO
C* z-center of the strip
C         ZSTR=ZSPEC+(3.5-NPLANE)*2.*MGPL(3)+(3-NPLANE)*2.*TRPB(3)+
C     +   2.*TRCP(3)+2.*TRRP(3)+2.*TSKP(3)
C+SELF
#if !defined(GPAMELA_PAOLO)
*         XSTR=XYZMO(1)
* Z-coord. of the implanted strip. The upper is the Y view
         ZSTR=XYZMO(3)-TSPA(3)
#endif
* magnetic field effect
         XI=XI-ABS(ZI-ZSTR)*0.2475*BB(2)/15.
         XO=XO-ABS(ZO-ZSTR)*0.2475*BB(2)/15.
*
C       PRINT *,ZSPEC,NPLANE,MGPL(3),TRPB(3),TRCP(3),TRRP(3),TSKP(3)
* charge definition for the strip IS
         CH_DIG(IS)=0.
* loop on the hit path division
         DO IH=1,NSEP
* x-coor. hit division
            XH=XI+(IH-0.5)*(XO-XI)/FLOAT(NSEP)
* Z-coor. hit division
            ZH=ZI+(IH-0.5)*(ZO-ZI)/FLOAT(NSEP)
            IX=1+INT(ABS(XH-XSTR)/STEPX)
            IZ=1+INT(ABS(ZH-ZSTR)/STEPZ)
C        PRINT *,'IX,IZ',IX,IZ,XH,XSTR,ZH,ZSTR
C        PRINT *,IS,IH,IX,IZ,CHARX(IX,IZ)
            IF(IX.LE.NPCHXY.AND.IZ.LE.NPCHZ) THEN
               IF(MOD(IS,2).EQ.1) THEN
                  CH_DIG(IS)=CH_DIG(IS)+DDE*CHARX(IX,IZ)
               ELSE
                  CH_DIG(IS)=CH_DIG(IS)+DDE*CHARX_I(IX,IZ)
               ENDIF
            ENDIF
         ENDDO
      END DO
* loop for C-coupling
      DO IS=NSTR_MIN,NSTR_MAX
* X-view has double strip (ISR is the read strip)
         ISR=INT((IS+1)/2)
         CHARGE=0.
         DO ISSR=ISR-2,ISR+2,1
            IF(2*ISSR-2.GE.NSTR_MIN.AND.2*ISSR.LE.NSTR_MAX) CHARGE=
     +      CHARGE+CCX(1+IABS(ISSR-ISR))* (CH_DIG(2*ISSR-1)+CH_DIG(2*
     +      ISSR)/2+CH_DIG(2*ISSR-2)/2)
         ENDDO
C       CHARGE=CH_DIG(IS)
*
* digitaz record
         NEWD=1
* X-view has double strip (ISR is the read strip)
         IF(LKDN.NE.0) THEN
            DO I=1,LKDN
               IF(NPLANE.EQ.LKD1(I).AND.NSIL.EQ.LKD2(I)
     +         .AND.ISR.EQ.LKD3(I)) THEN
                  NEWD=0
                  RLKD4(I)=RLKD4(I)+CHARGE
               ENDIF
            ENDDO
         ENDIF
         IF(NEWD.EQ.1) THEN
            LKDN=LKDN+1
            IF (LKDN.GT.NDMSPE) THEN
               WRITE(10000,CHMAIL) LKDN
               CALL GMAIL(1,0)
               GO TO 30
            ENDIF
            LKD1(LKDN)=NPLANE
            LKD2(LKDN)=NSIL
            LKD3(LKDN)=ISR
            RLKD4(LKDN)=CHARGE
            LKDH(LKDN)=NHIT
         ENDIF
      ENDDO
   10 CONTINUE
*
* Y-view digitalization
*
      NSEP=MAX0(NSEPY,NSEPZ)
      DDE=DE/NSEP
C+SELF,IF=PAOLO
C      NSIL=INT(((YI+YO)/2.+2.*TSPA(2))/(2.*TSPA(2)))+4
C      IF(NSIL.LT.4.OR.NSIL.GT.5) GOTO 20
C      YSIL_I=MOD(YI+2.*TSPA(2),2.*TSPA(2))
C      YSIL_O=MOD(YO+2.*TSPA(2),2.*TSPA(2))
C      YBORD=(2.*TSPA(2)-NYSTR*WYSTR)/2.
C      NSTR_I=INT((YSIL_I-YBORD)/WYSTR)+1
C      NSTR_O=INT((YSIL_O-YBORD)/WYSTR)+1
C      NSTR_MIN=MIN0(NSTR_I,NSTR_O)-5
C      IF(NSTR_MIN.LT.1) NSTR_MIN=1
C      NSTR_MAX=MAX0(NSTR_I,NSTR_O)+5
C      IF(NSTR_MAX.GT.NYSTR) NSTR_MAX=NYSTR
C+SELF.
#if !defined(GPAMELA_PAOLO)
      NSIL = 5 - (NUMV(NLEV)-1)
      IF(NSIL.LT.4.OR.NSIL.GT.5) GOTO 20
* Calculate the strip number now
      NSTR_I=INT((XYZISI(2))/WYSTR)+1
      NSTR_O=INT((XYZOSI(2))/WYSTR)+1
      NSTR_MIN=MIN0(NSTR_I,NSTR_O)-5
* minimum strip to be calculated
      IF(NSTR_MIN.LT.1) NSTR_MIN=1
      NSTR_MAX=MAX0(NSTR_I,NSTR_O)+5
* maximum strip to be calculated
      IF(NSTR_MAX.GT.NYSTR) NSTR_MAX=NYSTR
#endif
      DO IS=NSTR_MIN,NSTR_MAX
         YSTR=YBORD+(IS-0.5)*WYSTR+(NSIL-5)*2.*TSPA(2)
C+SELF,IF=PAOLO
C         ZSTR=ZSPEC+(3.5-NPLANE)*2.*MGPL(3)+(3-NPLANE)*2.*TRPB(3)+
C     +   2.*TRCP(3)+2.*TRRP(3)+2.*TSKP(3)+2.*TSPA(3)
C+SELF
#if !defined(GPAMELA_PAOLO)
*         YSTR=XYZMO(2)
* Z-coord. of the implanted strip. The lower is the X view
         ZSTR=XYZMO(3) + TSPA(3)
#endif
* magnetic field effect
         XI=XI+ABS(ZI-ZSTR)*0.04716*BB(2)/15.
         XO=XO+ABS(ZO-ZSTR)*0.04716*BB(2)/15.
*
         CH_DIG(IS)=0.
         DO IH=1,NSEP
            YH=YI+(IH-0.5)*(YO-YI)/NSEP
            ZH=ZI+(IH-0.5)*(ZO-ZI)/NSEP
            IY=1+INT(ABS(YH-YSTR)/STEPY)
            IZ=1+INT(ABS(ZH-ZSTR)/STEPZ)
C        PRINT *,IS,IH,IY,IZ,CHARY(IY,IZ)
            IF(IY.LE.NPCHXY.AND.IZ.LE.NPCHZ) CH_DIG(IS)=CH_DIG(IS)+DDE*
     +      CHARY(IY,IZ)
         ENDDO
      ENDDO
      DO IS=NSTR_MIN,NSTR_MAX
         CHARGE=0.
         DO ISS=IS-2,IS+2,1
            IF(ISS.GE.NSTR_MIN.AND.ISS.LE.NSTR_MAX) CHARGE=
     +      CHARGE+CCY(1+IABS(ISS-IS))*CH_DIG(ISS)
         ENDDO
*
* digitaz record
         NEWD=1
* the Y view paddles are connected
         IF(NSIL.EQ.5) NSIL=4
         IF(LKDN.NE.0) THEN
            DO I=1,LKDN
               IF(NPLANE.EQ.LKD1(I).AND.NSIL.EQ.LKD2(I)
     +         .AND.IS.EQ.LKD3(I)) THEN
                  NEWD=0
                  RLKD4(I)=RLKD4(I)+CHARGE
               ENDIF
            ENDDO
         ENDIF
         IF(NEWD.EQ.1) THEN
            LKDN=LKDN+1
            IF (LKDN.GT.NDMSPE) THEN
               WRITE(CHMAIL,10000) LKDN
               CALL GMAIL(1,0)
               GO TO 30
            ENDIF
            LKD1(LKDN)=NPLANE
            LKD2(LKDN)=NSIL
            LKD3(LKDN)=IS
            RLKD4(LKDN)=CHARGE
            LKDH(LKDN)=NHIT
         ENDIF
      ENDDO
   20 CONTINUE
*
10000 FORMAT(' GPSPEHD: too much digitalization= ',I5)
   30 RETURN
      END
