C   15/03/79 202262144  MEMBER NAME  LGCLBK   (S)           FORTRAN
C   13/03/79 C9031501   MEMBER NAME  LGCLPBN  (S)           FORTRAN
      SUBROUTINE LGCLBK
C     MODIFIED TO CALCULATE THE CLUSTER SHAPE
C     BY R.EICHLER       15-03-78  12:00
C     LAST MODIFICATION   16-08-79  19:30   Y.WATANABE
C     LAST MODIFICATION   26-02-82  21:45   S.YAMADA
C                                           POSITION OF COSD()=0.
C     ORIGINAL S.YAMADA 09-10-78  15:55 VERSION 4(L.M. 08-02-79,23:10)
C
C---- LEAD GLASS CLUSTERS IN THE BARREL PART ARE PROCESSED.
C---- TOTAL PULSE HEIGHT,AVERAGE POSITION AND SPATIAL SPREAD ARE
C     CALCULATED FOR EACH CLUSTER.
C
C---- MODIFICATION FOR VERSION 4;
C     ONLY THESE VARIABLES WHICH DO NOT NEED INNER TRACK INF.ARE FILLED.
C
      IMPLICIT INTEGER*2 (H)
C
      COMMON /CLGDMS/ X0,RADIUS(6),RADSX0(6),THX0(4),
     $                ZEND(2),ZENDX0(2),ZWID(2),ZGAP(2),PHWID(2),
     $                ZECAP(4),ZECAPX(4),THECPX(2)
      COMMON /CLGPRM/ ITHADC,MAXCLS
      COMMON/LGSHP/WW(50),ZZ(50),XX(50),EW(2),AVPH,AVZ,SIGPH,SIGZ,EV,KK
C------------------------------ ADDED BY S.K. --- 17-11-79 -------------
      COMMON /EIGVEC/EIGV11,EIGV12,EIGV21,EIGV22
      COMMON /EIGAB/COSD(2),COSA(2),COSB(2),
     &          EIGA1,EIGA2,EIGB1,EIGB2,ENECLA,ENECLB,NHITA,NHITB,DNORM
C     COMMON /CLVEC/V11,V12,V21,V22
C-----------------------------------------------------------------------
C
#include "clgwork1.for"
      EQUIVALENCE (NCLBRL,NCLBEC(1))
C
      COMMON /CLGMSB/ AMSGVL(5)
C-- MODIFIED BY S.KOMAMIYA -------------
      COMMON /MAXCLS/EIG1,EIG2,EIGSUM,EIGDIV,MAXHIT,ECMAX,MAXCL
      DIMENSION NHIT(20)
      DIMENSION EIGN(2,20)
      DIMENSION COS(2,20),PHAI(20),ZET(20)
C     DIMENSION EVC11(20)
C     DIMENSION EVC12(20)
C     DIMENSION EVC21(20)
C     DIMENSION EVC22(20)
C-------------------------------
C
C     ICLSPR(1,I)=JBC, 0 FOR BARREL, -1 FOR BOTTOM, 1 FOR TOP
C     CLSPRP(2,I)=ENERGY IN GEV
C     CLSPRP(3,I)=SIGMA(ENERGY)
C     CLSPRP(4,I)=WEIGHTED AVERAGE PHI
C     CLSPRP(5,I)=WEIGHTED AVERAGE Z
C     CLSPRP(6,I)=SIGMA PHI (WEIGHTED)
C     CLSPRP(7,I)=SIGMA Z (WEIGHTED)
C     ICLSPR(8,I)=CHARGE,IF THERE IS A CHARGED NEIGHBOUR
C                  TRACK,ITS CHARGE IS COPIED.
C     CLSPRP(9-11,I)=DIRECTION COSIGNS CORRECTED FOR SHOWER DEV.
C     CLSPRP(12-14,I)=CLUSTER SHAPE: ELLIPS EIGENVALUES AND DIRECTION
C     CLSPRP(15,I)=FRACTION OF ENERGY IN EDGE COUNTERS IN THE CLUSTER
C
C
      COMMON /CLGZPT/ YETOT,Y(32),WEIGHT(32),YFIT(32),DERIV(32),PHI(84)
C---- THIS COMMON IS USED TO GET THE (Z,PHI)-POSITION OF THE CLUSTER.
C     SEE LGAVRZ AND LGAVRP.
C
      DATA PHUNIT/0.07479983/
      COMMON /CLGVRN/ NVRSN(20)
      DATA NVCODE/478110108/
      NVRSN(9) = NVCODE
C
      ZSTEP=ZWID(1)+ZGAP(1)
      Z0=ZEND(1)+0.5*ZWID(1)
C
      IF(NCLBRL.LE.0) RETURN
      N2=NCLBRL
CC    IF(N2 .GT. MAXCLS) N2=MAXCLS
      IF(N2 .GT. 20 ) N2=20
        DO 1 N=1,N2
        ICLSPR(1,N) = 0
C
        NS = HMAPCL(1,N)
        NL = HMAPCL(2,N)
        WSUM = 0.
        IDO = HLGADC(1,NS)
        IPHO = IDO/32
C--------------------------------------------- ADDED BY S.K. -----------
        PHIO=FLOAT(IPHO)*PHUNIT*RADIUS(3)
C-----------------------------------------------------------------------
        NPMIN = 42
C----   CLEAR Y
        CALL SETSL(Y,0,848,0)
        KK=NL-NS+1
C---    CHECK NO. OF HITS IN THE CLUSTER
        IF(KK.LE.50) GO TO 3
        KK = 50
C
    3     DO 2 NN=NS,NL
          ID = HLGADC(1,NN)
          IPH = ID/32
          IZ = ID-32*IPH
C
C----     ANGLE PHI IS CONSIDERED WITH RESPECT TO THE FIRST BLOCK
C         TO AVOID ANGLE UNCERTAINTY.
          IPHD = IPH-IPHO
          IF(IPHD.GT.42) IPHD = IPHD-84
          IF(IPHD.LE.-42) IPHD = IPHD+84
          NPHD = IPHD+42
C----     CHECK THE MIN.PHI COUNTER #.
          IF(NPHD.LT.NPMIN) NPMIN = NPHD
C
C----     PULSE HEIGHT OF THE NN-TH HIT IS CALIBRATED IN MEV.
C
          W = HLGADC(2,NN)
C
C----     ADD W
          WSUM = WSUM+W
C
          IF(W) 90,9,9
C----     NEGATIVE ENERGY,   CLEARED TO PROTECT 'LGEIGN'.
   90     AMSGVL(1) = W
          CALL LGMESG( 6, 2)
          W = 0.
C
C
C----   FILL ARRAYS FOR CLUSTER SHAPE ANALYSIS
    9     NM=NN-NS+1
          IF(NM.GT.50) GO TO 4
          WW(NM)=W
          ZZ(NM)=IZ*ZSTEP+Z0
          XX(NM)=IPHD*PHUNIT*RADIUS(3)
C
C----     CHECK CORNER HIT
    4     CALL LGCRNR( 0, IZ, IPHI, ICLS)
          IF(ICLS.EQ.1) CLSPRP(15,N)=CLSPRP(15,N)+W
C
    2     CONTINUE
C
C
C----   CHANGE ENERGY UNIT TO GEV.
        WSUM = WSUM*0.001
C----   CALCULATE EIGENVALUES AND 3.MOMENTS OF CLUSTER
C-- MODIVIED BY S.KOMAMIYA -----
        CALL LGEIGK(1.)
C-- DELETED BY S.KOMAMIYA ------
C       CHANGE SIGZ AND SIGPH IN LINEAR DIMENSION.(Y.WATANABE 8/9/79)
C       IF(SIGZ.GT.0.) SIGZ=SQRT(SIGZ)
C       IF(ABS(AVZ-ZZ(1)).LT.0.1) SIGZ=ZWID(1)/3.4641
C       (IF ONLY ONE ROW OF THAT DIMENSION FIRE,
C       SIGMA=COUNTER WIDTH/2/SQRT(3).
C       IF(SIGPH.GT.0.) SIGPH=SQRT(SIGPH)/RADIUS(3)
C       IF(ABS(AVPH-XX(1)).LT.0.1) SIGPH=PHWID(1)/3.4641/RADIUS(3)
C
C----   CALCULATE SIG(ENERGY)
        CALL LGEERR(0,WSUM,SIGENG)
        CLSPRP(2,N) = WSUM
        CLSPRP(3,N) = SIGENG
        CLSPRP(4,N) = AVPH/RADIUS(3)+IPHO*PHUNIT+0.5*PHUNIT
        CLSPRP(5,N) = AVZ
        CLSPRP(6,N) = SIGPH
        CLSPRP(7,N) = SIGZ
C-------------------------------------------- MODIFIED BY S.K. ---------
        IF(EW(1).EQ.0.) GO TO 111
        CLSPRP(12,N)=EW(2)/EW(1)
        GO TO 112
111     CLSPRP(12,N)=0.
112     CLSPRP(13,N)=EW(1)+EW(2)
C-----------------------------------------------------------------------
        CLSPRP(14,N)=EV
        CLSPRP(15,N)=0.001*CLSPRP(15,N)/WSUM
C
        ETOT(2) = ETOT(2)+WSUM
C
C---------------------------------------------- MODIFIED BY S.K. -------
      NHIT(N) = NL-NS+1
      EIGN(1,N)=EW(1)
      EIGN(2,N)=EW(2)
C-- CALCULATE DIRECTION OF THE EIGEN VECTOR ( BELONGS TO EW(1) )
      EIGV1=SQRT(EIGV11**2+EIGV12**2)
      IF(EIGV1.LE.0.) GO TO 500
      COS(1,N)=EIGV11/EIGV1
      COS(2,N)=EIGV12/EIGV1
      GO TO 501
  500 COS(1,N)=0.
      COS(2,N)=0.
  501 CONTINUE
C
C-- CALCULATE PHI- AND Z-COORD. OF THE CLUSTER
C--------------------------- UNWEIGHTED SUM ----------------------------
      PHAI(N)=AVPH+PHIO
      ZET(N)=AVZ
C-----------------------------------------------------------------------
C     EVC11(N)=EIGV11
C     EVC12(N)=EIGV12
C     EVC21(N)=EIGV21
C     EVC22(N)=EIGV22
C-----------------------------------------------------------------------
    1 CONTINUE
C
      ECMAX=0.
      DO 10 L=1,N2
      IF(CLSPRP(2,L).LT.ECMAX) GO TO 10
      MAXCL=L
      ECMAX=CLSPRP(2,L)
   10 CONTINUE
      EIG1=EIGN(1,MAXCL)
      EIG2=EIGN(2,MAXCL)
      MAXHIT=NHIT(MAXCL)
      EIGDIV=CLSPRP(12,MAXCL)
      EIGSUM=CLSPRP(13,MAXCL)
C     V11=EVC11(MAXCL)
C     V12=EVC12(MAXCL)
C     V21=EVC21(MAXCL)
C     V22=EVC22(MAXCL)
C-----------------------------------------------------------------------
      NCLSA=0
      NCLSB=0
      ENECLA=0.
      ENECLB=0.
      DO 700 NN=1,N2
      IF(CLSPRP(2,NN).LE.ENECLA) GO TO 710
      ENECLB=ENECLA
      NCLSB=NCLSA
      ENECLA=CLSPRP(2,NN)
      NCLSA=NN
      GO TO 700
  710 IF(CLSPRP(2,NN).LE.ENECLB) GO TO 700
      ENECLB=CLSPRP(2,NN)
      NCLSB=NN
  700 CONTINUE
C-----------------------------------------------------------------------
      EIGA1=0.
      EIGA2=0.
      EIGB1=0.
      EIGB2=0.
      NHITA=0
      NHITB=0
      COSA(1)=0.
      COSA(2)=0.
      COSB(1)=0.
      COSB(2)=0.
      COSD(1)=0.
      COSD(2)=0.
      DNORM=0.
      IF(NCLSA.EQ.0) GO TO 800
      EIGA1=EIGN(1,NCLSA)
      EIGA2=EIGN(2,NCLSA)
      NHITA=NHIT(NCLSA)
      COSA(1)=COS(1,NCLSA)
      COSA(2)=COS(2,NCLSA)
      IF(NCLSB.EQ.0) GO TO 800
      EIGB1=EIGN(1,NCLSB)
      EIGB2=EIGN(2,NCLSB)
      NHITB=NHIT(NCLSB)
      COSB(1)=COS(1,NCLSB)
      COSB(2)=COS(2,NCLSB)
C---- CALCULATE DIRECTION COSINE OF THE ANGLE BTW TWO CLUSTERS ---------
      DPHI=PHAI(NCLSA)-PHAI(NCLSB)
      PHUNI1=PHUNIT*RADIUS(3)
      IF(DPHI/PHUNI1.GT.42.) DPHI=DPHI-84.*PHUNI1
      IF(DPHI/PHUNI1.LT.-42.) DPHI=DPHI+84.*PHUNI1
      DZET=ZET(NCLSA)-ZET(NCLSB)
      DNORM=SQRT(DPHI**2+DZET**2)
      IF(DNORM.LE.0.) GO TO 600
      COSD(1)=DPHI/DNORM
      COSD(2)=DZET/DNORM
  600 CONTINUE
C------------- TEST FOR REJECTING TWO RELATED COSMIC CLUSTERS----------
  800 CONTINUE
  300 CONTINUE
      RETURN
      END
