C   14/03/84            MEMBER NAME  GGCLPC   (S)           FORTRAN
C   29/09/78 006261425  MEMBER NAME  GGCLPC   (S)           FORTRAN
      SUBROUTINE GGCLPC
C
C     H.WRIEDT      29-11-78
C     LAST MODIFICATION 26-06-80  14:10
C
C     ORIGINALLY WRITTEN AS 'LGCLPC' BY S.YAMADA
C     AS TO BE FOUND ON YAMLGSRC.S
C
C
C---- LEAD GLASS CLUSTERS IN THE FORWARD DETECTORS ARE PROCESSED.
C---- TOTAL PULSE HEIGHT,AVERAGE POSITION AND SPATIAL SPREAD ARE
C     CALCULATED FOR EACH CLUSTER.
C
      IMPLICIT INTEGER*4 (G), INTEGER*2 (H)
      LOGICAL LSWTCH, LACOL
C
      COMMON /CGGDMS/ X0,ZECAP(4),ZECAPX(4),THECPX(2),
     1                ECXLST(46), ECYLST(46)
      COMMON /CGGPRM/ ITHADC,MAXCLS,CLBCNS
C
      COMMON /CALICO/ FAKTOR(191)
C
      COMMON /CMSGCT/ MSGDUM(20),GGMMAX,GGMSG1,GGMSG2
      COMMON /CWORK/ IWORK1(42),LNG,HPOINT(4),HGGADC(2,192),IWORK2(263),
     &               ID,HTRACK(3),HCLST,HCLBEC(2),HNCLST,HCOL,HTYPE,HER,
     &               HCORR,HPBLAT,HDUMIN(2),HSCAL(36),HWORD(2),ACOLAN,
     &               ETOT(3),ENTOT(3),IWORK3(16),HCLMAP(2,51),
     &               CLSPRP(13,51)
C
      COMMON /CPLADC/ HADCBL(128),HADCCL(128),HADCZM(128),HADCZP(128),
     &                HADCEV(128),ACOLPL(102)
C
      DIMENSION HCLSPR(26,51),COSZM(3),COSZP(3)
      EQUIVALENCE (HCLSPR(1,1),CLSPRP(1,1))
C
C     HCLSPR(1,I)=NO. OF THE CLUSTER
C     HCLSPR(2,I)=JBC, +1 FOR +Z-DIRECTION (TOWARDS MARK J)
C                      -1 FOR -Z-DIRECTION (TOWARDS PLUTO/CELLO)
C                       0 FOR ILLEGAL COMPONENT NO.
C     HCLSPR(3,I)=NUMBER OF TRACKS WHICH GO THROUGH THIS CLUSTER
C     HCLSPR(4,I)=NO. OF THE FIRST TRACK WHICH GOES THROUGH THIS
C                 CLUSTER
C     HCLSPR(5,I)=NO. OF THE COLLINEAR CLUSTER
C     CLSPRP(5,I)=CLUSTER SHOWER ENERGY IN GEV
C     CLSPRP(6,I)=EXPECTED ERROR OF THE ENERGY
C     CLSPRP(7,I)=WEIGHTED AVERAGE X
C     CLSPRP(8,I)=WEIGHTED AVERAGE Y
C     CLSPRP(9,I)=SIGMA X (WEIGHTED)
C     CLSPRP(10,I)=SIGMA Y (WEIGHTED)
C     CLSPRP(11-13,I)=DIRECTION COSINES OF THE CLUSTER
C
      COMMON/ CGGVRN/ NVRSN(20)
      DATA NVCODE/278113010/
      NVRSN(10) = NVCODE
C
   10 IF (HCLST.LE.0) RETURN
C---  INITIALIZE ADC-SUMS
      ISUMZM = 0
      ISUMZP = 0
C
        DO 1 N=1,HCLST
C
        NS = HCLMAP(1,N)
        NL = HCLMAP(2,N)
        WSUM = 0.
        WTSM = 0.
        XSUM = 0.
        YSUM = 0.
        X2SM = 0.
        Y2SUM = 0.
C---  INITIALIZE ADC-SUM
        ISUMCL = 0
C
          DO 2 NN=NS,NL
          ID = HGGADC(1,NN)
          CALIF = FAKTOR(ID)
          IPART = -1
          IF (ID.EQ.47.OR.ID.EQ.48.OR.ID.EQ.95.OR.ID.EQ.96) GOTO 89
          IF(ID.LE.96) GO TO 3
          ID = ID-96
          IPART = +1
          IF (ID.EQ.47.OR.ID.EQ.48.OR.ID.GE.95) GOTO 89
    3     MB = (ID-1)/48
          IF (MB.EQ.1) ID = 95 - ID
          U = ECXLST(ID)
          V = ECYLST(ID)
          MB = (1-2*MB)*IPART
          IF (MB) 5,89,4
    4     X = U
          Y = -V
          GO TO 8
    5     X = -U
          Y = V
    8     W = HGGADC(2,NN)*CLBCNS*CALIF
          HGGADC(2,NN) = W*1000
C
C---  FILL ADC-PLOT FOR BLOCK-SUM
C     ASSUMED: 12-BIT ADCS
          IDAT = (HGGADC(2,NN)+31)/32
          HADCBL(IDAT) = HADCBL(IDAT)+1
          ISUMCL = ISUMCL + IDAT
C
C---- WEIGHT IS CHANGED TO SQRT(W)
          WT = SQRT(W)
C----     ADD W,X,X2,Y,Z2
          WSUM = WSUM+W
          WTSM = WTSM+WT
          WX = WT*X
          XSUM = XSUM+WX
          X2SM = X2SM+WX*X
          YW = Y*WT
          YSUM = YSUM+YW
          Y2SUM = Y2SUM+YW*Y
          GOTO 2
C-----    ILLEGAL COMPONENT NO.
 89       CONTINUE
          IF(GGMSG2) 12,90,90
   90     GGMSG2 = GGMSG2+1
          IF(GGMSG2.LE.GGMMAX) WRITE(6,600) NN,ID,HGGADC(1,NN)
  600     FORMAT(' LG FORWARD DETECTOR FORMAT ERROR,  ID,HGGADC(1,',I4,
     &           ')=',2I5)
   12     IF (HGGADC(2,NN).GT.10) IPART = 0
          IF (HGGADC(2,NN).LE.10) CALL GGRESQ(NN,N,&10)
C-----    BAD DATA IS SKIPPED.
    2    CONTINUE
        IF (IPART.NE.0) GOTO 6
C---  FILL CLSPRP FOR CLUSTERS FROM ILLEGAL COMPONENTS
        HCLSPR(1,N) = N
        HCLSPR(2,N) = IPART
        CLSPRP(5,N) = 0.
        CLSPRP(7,N) = 0.
        CLSPRP(8,N) = 0.
        CLSPRP(9,N) = 0.
        CLSPRP(10,N) = 0.
        CLSPRP(11,N) = 0.
        CLSPRP(12,N) = 0.
        CLSPRP(13,N) = 0.
        GOTO 1
C---  FILL ADC-PLOT FOR CLUSTER-SUM
    6   IF (ISUMCL.GT.128) ISUMCL=128
        HADCCL(ISUMCL) = HADCCL(ISUMCL)+1
        IF (IPART.EQ.-1) ISUMZM = ISUMZM+ISUMCL
        IF (IPART.EQ.+1) ISUMZP = ISUMZP+ISUMCL
C
        IF (WTSM.GT.0.) GOTO 13
        WRITE(6,610) ((HGGADC(K,J),K=1,2),J=NS,NL)
        WRITE(6,611) WTSM,XSUM,YSUM
  610   FORMAT(/' *** ENERGY-SUM = 0 ***'/20(2X,I4))
  611   FORMAT(3(2X,I4))
        AVX = 0.
        SIGX = 0.
        AVY = 0.
        SIGY = 0.
        GOTO 14
   13   AVX = XSUM/WTSM
        SIGX = X2SM/WTSM-AVX*AVX
        AVY = YSUM/WTSM
        SIGY = Y2SUM/WTSM-AVY*AVY
C---  SUM OF TOTAL ENERGIES IN EACH DETECTOR HALF
   14   XIND = IPART/2. + 2.6
        INDEX = XIND
        ETOT(INDEX) = ETOT(INDEX) + WSUM
C---  FILL CLSPRP
        HCLSPR(1,N) = N
        HCLSPR(2,N) = IPART
        CLSPRP(5,N) = WSUM
        CLSPRP(7,N) = AVX
        CLSPRP(8,N) = AVY
        CLSPRP(9,N) = SIGX
        CLSPRP(10,N) = SIGY
C       CALCULATE THETA AND PHI
        ZVALUE = ZECAP(IPART+2)
        TR = SQRT(AVX*AVX+AVY*AVY+ZVALUE*ZVALUE)
        CLSPRP(11,N) = AVX/TR
        CLSPRP(12,N) = AVY/TR
        CLSPRP(13,N) = ZVALUE/TR
C---  SEARCH CLUSTER WITH MAXIMUM ENERGY IN EACH DETECTOR HALF
        IF (N.GE.2) GOTO 20
        LSWTCH = .FALSE.
        LACOL = .FALSE.
        IF (IPART.EQ.1) GOTO 1
C---  -Z-PART
        LACOL = .TRUE.
          DO 15 I = 1,3
   15     COSZM(I) = CLSPRP(I+10,1)
        ENMAX = CLSPRP(5,1)
        GOTO 1
C
   20   IF (.NOT.LACOL) GOTO 1
        IF (LSWTCH) GOTO 40
        IF (IPART.EQ.1) GOTO 30
        IF (CLSPRP(5,N).LE.ENMAX) GOTO 1
          DO 25 I = 1,3
   25     COSZM(I) = CLSPRP(I+10,N)
        ENMAX = CLSPRP(5,N)
        GOTO 1
C---  +Z-PART
   30   LSWTCH = .TRUE.
          DO 35 I = 1,3
   35     COSZP(I) = CLSPRP(I+10,N)
        INDEX = N
        CALL GGCOLS(INDEX)
        ENMAX = CLSPRP(5,N)
        GOTO 1
C
   40   IF (CLSPRP(5,N).LE.ENMAX) GOTO 46
          DO 45 I = 1,3
   45     COSZP(I) = CLSPRP(I+10,N)
        ENMAX = CLSPRP(5,N)
   46   INDEX = N
        CALL GGCOLS(INDEX)
    1   CONTINUE
C
C---  ADD TOTAL ENERGIES OF BOTH DETECTOR HALVES
      ETOT(1) = ETOT(2) + ETOT(3)
C---  ENERGY OF CLUSTERS WITHOUT TRACKS OR LUMONITORS
        DO 50 I = 1,3
   50   ENTOT(I) = ETOT(I)
C---  FILL ADC-PLOT FOR DETECTOR HALVES-SUM AND EVENT-SUM
      ISUMEV = ISUMZM + ISUMZP
      IF (ISUMZM.GT.128) ISUMZM=128
      IF (ISUMZM.GT.0) HADCZM(ISUMZM) = HADCZM(ISUMZM) + 1
      IF (ISUMZP.GT.128) ISUMZP=128
      IF (ISUMZP.GT.0) HADCZP(ISUMZP) = HADCZP(ISUMZP) + 1
      IF (ISUMEV.GT.128) ISUMEV=128
      HADCEV(ISUMEV) = HADCEV(ISUMEV) + 1
C
C---  CALCULATE ACOLLINEARITY ANGLE
      IF (LACOL.AND.LSWTCH) ACOLAN = ACOLL(COSZM,COSZP)
C     ACOLLINEARITY ANGLE DISTRIBUTION
      CALL SORT(0.,ACOLAN,ACOLPL,0.002,100,1)
  100 RETURN
      END
      FUNCTION ACOLL(COSZM,COSZP)
C---  CALCULATE ACOLLINEARITY ANGLE BETWEEN TWO ANGLES, GIVEN BY THEIR
C     DIRECTION COSINES, COSZM AND COSZP
C
C     H.WRIEDT      31.01.79         18:40
C     LAST MODIFICATION    01.02.79     19:05
C
      DIMENSION COSZM(3), COSZP(3)
C
      SUM = 0
        DO 10 I = 1,3
  10    SUM = SUM - COSZM(I)*COSZP(I)
C
C***  CHECK COSINE RANGE
      IF (SUM.GE.-1. .AND. SUM.LE.1.) GOTO 20
      IF (SUM.GE.-1.01 .OR. SUM.LE.1.01) GOTO 15
      WRITE(6,100) SUM
 100  FORMAT(' ***** ERROR IN ACOLAN, COSINE OUT OF RANGE: ',E16.8)
      ACOLL = 9.99
      RETURN
C---  CURE ROUNDING ERROR
  15  SUM = SIGN(1.,SUM)
C
  20  ACOLL = ARCOS(SUM)
      RETURN
      END
      SUBROUTINE GGCOLS(INDEX)
C---  FIND COLLINEAR PAIRS OF CLUSTERS
C---  A PAIR OF CLUSTERS IS DEFINED AS "COLLINEAR" IF THE ACOLLINEARITY
C---  ANGLE BETWEEN THE TWO CLUSTERS IS LESS THAN 17.5 MRAD (ABOUT
C---  1 DEGREE)
C
C     H.WRIEDT     22.03.79       20:50
C     LAST MODIFICATION     11.04.79     11:10
C
      IMPLICIT INTEGER*2 (H)
C
      COMMON /CWORK/ IWORK(503),HCLBEC(2),HNCLST,HCOL,IWORK2(96),
     &               CLSPRP(13,51)
      DIMENSION HCLSPR(26,51),COSZM(3),COSZP(3)
      EQUIVALENCE (HCLSPR(1,1),CLSPRP(1,1))
C
        DO 5 J = 1,3
    5   COSZP(J) = CLSPRP(J+10,INDEX)
C
      N = HCLBEC(1)
        DO 10 I = 1,N
          DO 1 J = 1,3
    1     COSZM(J) = CLSPRP(J+10,I)
        IF (ACOLL(COSZM,COSZP).GT.17.5E-3) GOTO 10
        HCOL = HCOL + 1
        HCLSPR(5,INDEX) = HCLSPR(1,I)
        HCLSPR(5,I) = HCLSPR(1,INDEX)
        RETURN
   10   CONTINUE
      RETURN
      END
