***********************************************************************
*
* Subroutine READVME
* Routine to read a single event: use C subroutine readevnt.
* Use conditional compilation with cpp preprocessor, so this file must 
* end in .F and cannot be included with INCLUDE directive (use #include)
* Input:  IEV event # to read 
* Output: fill ADC and DSP variables
*
* Rewritten from a previous routine from an HP machine
* S. Piccardi Nov. 1998
*
* $Id: readvme.F,v 1.1 1999/12/15 10:05:56 piccardi Exp $
*
*****************************************************************************
      PROGRAM NNPLOT_NT5
c
      IMPLICIT NONE
C
      INCLUDE 'INTEST.TXT'
C
      CHARACTER*6 file_name, file_name2
      CHARACTER*40 file_name3
C      
      integer icount,icount2
c
      integer RUNERROR, PEDERROR
C     
C Normal variables definition
C
      INTEGER FFD
C
      INTEGER l, j, kk, m, mm, ll, il, ii, nn, ival, iev
      INTEGER i, istat, icycle, ierr, icheck, icheck2, icheck3 
      INTEGER LVMIN

      INTEGER*2 length

      integer IP

      INTEGER ipre, ipre2, jc
      INTEGER ic, k, jj
      INTEGER status
      INTEGER inf, isup, idet

      INTEGER LEVEL

      REAL ADCERR

C$$$      INTEGER*4 lenght(4)

      INTEGER lundata, iosop

      INTEGER*2 st, st2

      CHARACTER*2 nr
c
      integer*4 evnum, systime, len
      integer nev
C
      integer nset,nword(10) 
C
      INTEGER*2 STRINGA(200000) 
C
      REAL auto(11)

      integer lun(4)

      REAL SUM(6), SUM2(6), COUNT(6), 
     &     BAS(4,11,6), DIS(6), DDIS(6), BASE(2)
      REAL SSUM, SSUM2, CCOUNT, MME(2), SSIG

      REAL VAL

      REAL PED(4,11,96), GOOD(4,11,96)
      COMMON / PEDE / PED, GOOD

      REAL MIP1(11,96), MIP2(11,96), MIP3(11,96), MIP4(11,96)
      REAL SPOS1(3,11), SPOS2(3,11), SPOS3(3,11), SPOS4(3,11), DET

      REAL D1(11,96), D2(11,96), D3(11,96), D4(11,96)

      REAL DI, DS
      REAL ME(6), SIG(6)

      REAL TG(2)
      REAL SHIFT
      COMMON /SHIFT/ SHIFT

      REAL BAR(2,NPLA)
      INTEGER IBAR(2,NPLA)
      COMMON/ANGOLO/BAR,IBAR

      REAL DISTX, DISTY, Y(NPLA), YY(NPLA)

      REAL CX, CY
      COMMON/WHERE/CX,CY

      REAL RIG, PLANEMAX, RMASS
      COMMON/GENERAL/RIG,RMASS
  
      INTEGER IPLANE, NNX, NNY, INFX, INFY, ISUPX, ISUPY
      INTEGER NLK

      REAL RNSS, QTOTT, RQT

      real dedx1(11,96),dedx2(11,96),dedx3(11,96),dedx4(11,96)
      real qpl(96), ene(96), qst(96)

      INTEGER nin

      REAL CHECK
      COMMON / CH / CHECK

      REAL QMAXX, QMAXY
      REAL VMAX

      real nstrip, qtot, ncore, qcore, impx, impy, tanx, tany, nint
      REAL ncyl, qcyl, qtrack, qmax, cher, lead, nx22, qx22, qq(4)
      INTEGER event
      COMMON / Planes / nstrip, qtot, ncore, qcore, impx, impy, 
     &     tanx, tany, nint, ncyl, qcyl, qtrack, qmax, cher, lead, nx22, 
     &     qx22, qq, event

      INTEGER lrec
      PARAMETER (lrec=6144)

      REAL hmemor(3000000)
      COMMON /pawc/hmemor
      CALL HLIMIT(3000000)

C         
C Begin !
C
      RUNERROR = 0                ! error variable
      PEDERROR = 0

      PRINT *,'Data file to read (only data and number included)? '
      READ(*,905)file_name

      CALL PEDESTAL(PEDERROR)
      IF (PEDERROR.EQ.1) THEN
         PRINT *,'Problems with pedestal file'
         GOTO 50
      ENDIF
C
905   FORMAT(A6)

      lundata=33

      OPEN(UNIT=lundata,FILE='/wizard3/pamela/tb-0903/data/output_0309'
     +  //file_name//'.dat',STATUS='OLD',
     +  FORM='UNFORMATTED',ERR=50,IOSTAT=IOSOP)
C
C
C
      PRINT *,'File to save? '
      READ(*,904)file_name3

 904  FORMAT(A40)
C
C   Histos creation
C
      CALL HROPEN(59,'Event','/wizard3/pamela/tb-0903/paw/'//
     +                  file_name3,'n',lrec,istat)
      CALL HBNT(1,'Pamela Calo',' ')
      CALL HBSET('BSIZE',lrec,ierr)

***     /* Book ntuple variables */
      CALL HBNAME(1,'calo',nstrip,'nstrip:R,qtot:R,
     &     ncore:R,qcore:R,impx:R,impy:R,tanx:R,tany:R,
     &     nint:R,ncyl:R,qcyl:R,qtrack:R,qmax:R,
     &     cher:R,lead:R,nx22:R,qx22:R,qq(4):R,event:I')
C
 906  FORMAT(A2)
C
      OPEN(71,FILE='/wizard3/pamela/tb-0903/source/mip1.dat',
     &     STATUS='OLD',FORM='FORMATTED',ERR=50,IOSTAT=IOSOP)
      OPEN(72,FILE='/wizard3/pamela/tb-0903/source/mip2.dat',
     &     STATUS='OLD',FORM='FORMATTED',ERR=50,IOSTAT=IOSOP)
      OPEN(73,FILE='/wizard3/pamela/tb-0903/source/nmip3.dat',
     &     STATUS='OLD',FORM='FORMATTED',ERR=50,IOSTAT=IOSOP)
      OPEN(74,FILE='/wizard3/pamela/tb-0903/source/mip4.dat',
     &     STATUS='OLD',FORM='FORMATTED',ERR=50,IOSTAT=IOSOP)
C
      OPEN(81,FILE='/wizard3/pamela/tb-0903/source/spos1.dat',
     &     STATUS='OLD',FORM='FORMATTED',ERR=50,IOSTAT=IOSOP)
      OPEN(82,FILE='/wizard3/pamela/tb-0903/source/spos4.dat',
     &     STATUS='OLD',FORM='FORMATTED',ERR=50,IOSTAT=IOSOP)
      OPEN(83,FILE='/wizard3/pamela/tb-0903/source/spos3.dat',
     &     STATUS='OLD',FORM='FORMATTED',ERR=50,IOSTAT=IOSOP)
      OPEN(84,FILE='/wizard3/pamela/tb-0903/source/spos2.dat',
     &     STATUS='OLD',FORM='FORMATTED',ERR=50,IOSTAT=IOSOP)
C
      DO I = 1,11
         READ (71,102)(mip1(I,J),J=1,96)
         READ (72,102)(mip2(I,J),J=1,96)
         READ (73,102)(mip3(I,J),J=1,96)
         READ (74,102)(mip4(I,J),J=1,96)
      ENDDO
      DO I = 1,3
         READ (81,103)(SPOS1(I,J),J=1,11)
         READ (82,103)(SPOS2(I,J),J=1,11)
         READ (83,103)(SPOS3(I,J),J=1,11)
         READ (84,103)(SPOS4(I,J),J=1,11)
         print *,'I :',I,' J : 2 spos1 :',SPOS1(I,2)
      ENDDO
      CLOSE(21)
      CLOSE(71)
      CLOSE(72)
      CLOSE(73)
      CLOSE(74)
      CLOSE(81)
      CLOSE(82)
      CLOSE(83)
      CLOSE(84)

 102  FORMAT (4(2X,F9.3))
 103  FORMAT (11(2X,F7.5))
C
C Take the file descriptor number (unix)       
C
      PRINT *,'Momentum:'
      READ *,RIG
C
      PRINT *,'Mass of the particle:'
      READ *,RMASS
C
      PRINT *,'How many events to read:'
      READ *,nev
C
      ffd=FNum(lundata)          ! take the file descriptor number
C

C Read event: use C readevent routine
C
      icount = 0
      icount2 = 0
      event = 0
      iev = 0.
      CALL VZERO(BAS,66*4)
c
 10   continue
C
      event = event + 1
      ival = 0
      DO I = 1,11
         DO J = 1,96
            DEDX1(I,J) = 0.
            DEDX2(I,J) = 0.
            DEDX3(I,J) = 0.
            DEDX4(I,J) = 0.
         ENDDO
      ENDDO
C
      call readev(evnum,systime,len,stringa,RUNERROR,ffd)  ! XXX

      if (runerror.eq.-1.or.runerror.eq.1) goto 50
C
c      write (*,31) evnum
C      write (*,32) systime
C      write (*,33) len
C      print *,'event number:',evnum
C      print *,'systime:',systime
C      print *,'event length:',len
 31   format(1x,'event number :',1x,Z8)
 32   format(1x,'systime :',1x,Z8)
 33   format(1x,'event length :',1x,Z8)

C$$$      write(6,*)nset,nword
 
C$$      print *,'lenght 1 :',lenght(1)
C$$      print *,'lenght 2 :',lenght(2)
C$$      print *,'lenght 3 :',lenght(3)

      if (runerror.eq.-1.or.runerror.eq.1) goto 50

      ic = 0
 
      ic = ic + 4
      length = ic
      do kk = 1,4
C         if (kk.eq.3) goto 45
         ic = length + 1
C$$$         do k = ic-5,ic+5
         k = ic
C         write (*,41) stringa(k)
         st2 = IAND(stringa(k),'FF00'x)
         status = ISHFT(st2,-8)
C         write (*,41) status
         IF (KK.eq.1.and.status.ne.234) then
            print *,'Problems with view:',KK
            GOTO 45
         endif
         IF (KK.eq.2.and.status.ne.246) then
            print *,'Problems with view:',KK
            GOTO 45
         endif
         IF (KK.eq.3.and.status.ne.241) then
            print *,'Problems with view:',KK
            GOTO 45
         endif
         IF (KK.eq.4.and.status.ne.237) then
            print *,'Problems with view:',KK
            GOTO 45
         endif

C$$$         enddo
         if (status.ne.-30460) then
c            print *,'problems with DSP, event:',event
         endif
C
 41      format(1x,'status :',1x,Z4)
C
         ic = ic + 1
         lun(kk) = stringa(ic)
         length = length + (stringa(ic) + 2)
C         print *,'Lenght1 :',length
c         write (*,41) stringa(ic)
C         print *,'Lenght2 :',stringa(ic)

C$$$         if (stringa(ic).ne.1064) then
C$$$            print *,'problems with DSP'
C$$$         endif
C
         do i = 1,11
            ic = ic + 1
            auto(i) = stringa(ic)
         enddo
C
         do i = 1,11
            do j = 1,96
               ic = ic + 1
               if (kk.eq.1) then 
                  d1(i,j) = stringa(ic) - 
     &                 ped(kk,i,j)
                     if (i.eq.4.and.j.gt.64) then
c                        print *,'d1 :',d1(i,j),' ped :',ped(kk,i,j),
c     &                       ' stringa :',stringa(ic),' j :',j
                     endif
               endif
               if (kk.eq.2) then 
                  d2(i,j) = stringa(ic) - 
     &                 ped(kk,i,j)
               endif
               if (kk.eq.3) then 
                  d3(i,j) = stringa(ic) - 
     &                 ped(kk,i,j)
               endif
               if (kk.eq.4) then 
                  d4(i,j) = stringa(ic) - 
     &             ped(kk,i,j)
C$                  if (i.eq.2.and.j.gt.32.and.j.lt.49.and.event.eq.42) 
C$     &                 then
C$                        print *,'d4 :',d4(i,j),' ped :',ped(kk,i,j),
C$     &                       ' stringa :',stringa(ic),' j :',j
C$                     endif
               endif
            enddo
         enddo
C
         DO I = 1,11
C
            DO M = 1,6
               INF = (M - 1) * 16 + 1
               ISUP = M * 16
               CHECK = 0.
C$               IF (M.EQ.3.AND.I.EQ.2.AND.KK.EQ.4.and.event.eq.42) 
C$     &              CHECK = 1.
               IF (KK.EQ.1) CALL CALCBASE(KK,INF,ISUP,I,D1,BAS(KK,I,M))
c               IF (KK.EQ.1.and.i.eq.4) PRINT *,'BASE 1:',BAS(KK,I,M),
c     &              ' m :',m
               IF (KK.EQ.2) CALL CALCBASE(KK,INF,ISUP,I,D2,BAS(KK,I,M))
               IF (KK.EQ.3) CALL CALCBASE(KK,INF,ISUP,I,D3,BAS(KK,I,M))
c               IF (M.EQ.4.AND.I.EQ.4.AND.KK.EQ.4) 
c     &              CHECK = 1.
               IF (KK.EQ.4) CALL CALCBASE(KK,INF,ISUP,I,D4,BAS(KK,I,M))
            ENDDO
C     
C
            do j = 1,96
               IP = INT((J - 1) / 16) + 1
               IF (KK.EQ.1) THEN
                  IF (BAS(1,I,IP).GT.0..AND.MIP1(I,J).GT.0.) THEN
                     DEDX1(I,J) = (D1(I,J) - BAS(1,I,IP)) / MIP1(I,J)
C                     if (i.eq.1) print *,'dedx :',dedx1(I,j),
C     &                    ' base :',bas(1,i,ip),' mip :',mip1(i,j)
                     if (i.eq.4.and.j.gt.64) then
c                        print *,'dedx1 :',dedx1(i,j),' i :',i,' j :',j
                     endif
                  ELSE
                     DEDX1(I,J) = 0.
                  ENDIF
               ENDIF
               IF (KK.EQ.2) THEN
                  IF (BAS(2,I,IP).GT.0..AND.MIP2(I,J).GT.0.) THEN 
                     DEDX2(I,J) = (D2(I,J) - BAS(2,I,IP)) / MIP2(I,J) 
                  ELSE
                     DEDX2(I,J) = 0.
                  ENDIF
               ENDIF
c               IF (KK.EQ.3) THEN
c                  IF (BAS(3,I,IP).GT.0.) THEN
c                     DEDX3(I,J) = (D3(I,J) - BAS(3,I,IP)) / MIP3(I,J)
c                  ELSE
c                     DEDX3(I,J) = 0.
c                  ENDIF
c               ENDIF
               IF (KK.EQ.4) THEN
                  IF (BAS(4,I,IP).GT.0..AND.MIP4(I,J).GT.0.) THEN
                     DEDX4(I,J) = (D4(I,J) - BAS(4,I,IP)) / MIP4(I,J)
c                     IF (I.EQ.4.and.j.gt.48.and.j.lt.65) 
c                          PRINT *,'DEDX4 :',DEDX4(I,J),' BASE :',
c     &                    BAS(4,I,IP),' D4 :',D4(I,J),' J :',J
                  ELSE
                     DEDX4(I,J) = 0.
                  ENDIF
               ENDIF
            enddo
c     
         enddo
C
         ic = ic + 1
C
 45      continue
      enddo
C
      DO M = 1,4
         DO K = 1,11
C
            DO I = 1,96
               IF (M.EQ.1.AND.MIP1(K,I).NE.0.) QPL(I) = DEDX1(K,I)
               IF (M.EQ.2.AND.MIP2(K,I).NE.0.) QPL(I) = DEDX2(K,I)
               IF (M.EQ.3.AND.MIP3(K,I).NE.0.) QPL(I) = DEDX3(K,I)
               IF (M.EQ.4.AND.MIP4(K,I).NE.0.) QPL(I) = DEDX4(K,I)
               ENE(I) = QPL(I)
               QST(I) = QPL(I)
            ENDDO   
C
            DO I = 1,96
               IPRE = INT((I-1)/16) + 1
               JC = MOD(I,32)
               IF (QPL(I).GT.EMIN) THEN
C
                  INF = (IPRE - 1) * 16 + 1
                  ISUP = IPRE * 16
                  DO J = INF,ISUP
C$$$                     ENE(J) = ENE(J) + QPL(I) * 0.00963
                     ENE(J) = ENE(J) + QPL(I) * 0.00478
                  ENDDO
C
                  IF (MOD(IPRE,2).EQ.0) IPRE2 = IPRE - 1
                  IF (MOD(IPRE,2).EQ.1) IPRE2 = IPRE + 1
                  INF = (IPRE2 - 1) * 16 + 1
                  ISUP = IPRE2 * 16
C$$$                  DO J = INF,ISUP
C$$$                     ENE(J) = ENE(J) + QPL(I) * 0.00485
C$$$                  ENDDO
C
                  IF (JC.NE.1) ENE(I-1) = ENE(I-1) - QPL(I)*0.01581
                  IF (JC.NE.0) ENE(I+1) = ENE(I+1) - QPL(I)*0.01581
C     
               ENDIF
            ENDDO
C
            DO I = 1,96
               IPRE = INT((I-1)/16) + 1
               JC = MOD(I,32)
               IF (ENE(I).GT.EMIN) THEN
C     
                  INF = (IPRE - 1) * 16 + 1
                  ISUP = IPRE * 16
                  IDET = INT((I - 1) / 32) + 1
                  IF (M.EQ.1) DET = SPOS1(IDET,K)
                  IF (M.EQ.2) DET = SPOS2(IDET,K)
                  IF (M.EQ.3) DET = SPOS3(IDET,K)
                  IF (M.EQ.4) DET = SPOS4(IDET,K)
                  DO J = INF,ISUP
C$$$                     QST(J) = QST(J) + ENE(I) * 0.00963
                     QST(J) = QST(J) + ENE(I) * (0.00478 + 
     &                    DET)
                  ENDDO
C     
                  IF (MOD(IPRE,2).EQ.0) IPRE2 = IPRE - 1
                  IF (MOD(IPRE,2).EQ.1) IPRE2 = IPRE + 1
                  INF = (IPRE2 - 1) * 16 + 1
                  ISUP = IPRE2 * 16
                  DO J = INF,ISUP
                     QST(J) = QST(J) + ENE(I) * DET
                  ENDDO
C     
                  IF (JC.NE.1) QST(I-1) = QST(I-1) - ENE(I)*0.01581
                  IF (JC.NE.0) QST(I+1) = QST(I+1) - ENE(I)*0.01581
C     
               ENDIF
            ENDDO
C
            DO I = 1,96
C               IF (M.EQ.1.AND.MIP1(K,I).NE.0.) DEDX1(K,I) = QST(I)
c               IF (M.EQ.2.AND.MIP2(K,I).NE.0.) DEDX2(K,I) = QST(I)
c               IF (M.EQ.3.AND.MIP3(K,I).NE.0.) DEDX3(K,I) = QST(I)
c               IF (M.EQ.4.AND.MIP4(K,I).NE.0.) DEDX4(K,I) = QST(I)
               IF (M.EQ.4.AND.K.EQ.5) DEDX4(K,I) = 0.
            ENDDO
C
         ENDDO
C         
      ENDDO
C
      CALL VZERO(DEXY,2*LENSEV)
      CALL VZERO(QQ,4)
      NSTRIP = 0.
      QTOT = 0.
      NX22 = 0.
      QX22 = 0.
      DO I = 1,11
         DO J = 1,96
            ICHECK3 = 1
            IF (DEDX1(I,J).GT.EMIN)
     &           THEN
               DEXY(2,2*I-1,97-J) = DEDX1(I,J)
               NSTRIP = NSTRIP + 1.
               QTOT = QTOT + DEDX1(I,J)
               IF (I.NE.1) QQ(1) = QQ(1) + DEDX1(I,J)
c               if (i.eq.1) print *,'qtot1 :',qtot,' dedx :',dedx1(i,j)
            ENDIF
            IF (DEDX2(I,J).GT.EMIN) THEN
               DEXY(1,2*I-1,J) = DEDX2(I,J)
               NSTRIP = NSTRIP + 1
               QTOT = QTOT + DEDX2(I,J)
               IF (I.EQ.11) THEN
                  NX22 = NX22 + 1.
                  QX22 = QX22 + DEDX2(I,J)
               ENDIF
               if (i.lt.11) QQ(2) = QQ(2) + DEDX2(I,J)
c               print *,'qtot2 :',qtot
            ENDIF
            IF (DEDX3(I,J).GT.EMIN) THEN
C               DEXY(2,2*I,J) = DEDX3(I,J)
C               NSTRIP = NSTRIP + 1
C               QTOT = QTOT + DEDX3(I,J)
C               QQ(3) = QQ(3) + DEDX3(I,J)
c               print *,'qtot3 :',qtot
            ENDIF
            IF (DEDX4(I,J).GT.EMIN) THEN
               DEXY(1,2*I,J) = DEDX4(I,J)
               NSTRIP = NSTRIP + 1
               QTOT = QTOT + DEDX4(I,J)
               IF (I.LT.11) QQ(4) = QQ(4) + DEDX4(I,J)
c               if (i.eq.5) PRINT *,'DEDX4 :',DEDX4(I,J),' qtot :',
c     &              qtot,' J :',J,' I :',i
c               print *,'qtot4 :',qtot
            ENDIF
         ENDDO
      ENDDO
C
      CALL CLUSTER
C
      CALL DIRECTION(TG)
      IMPX = CX - 241. / 2.
      IMPY = CY - 241. / 2.
      SHIFT = -0.5
      CALL LASTRISCIA(CX,II)
      IMPX = II
      SHIFT = +0.5
      CALL LASTRISCIA(CY,II)
      IMPY = II
      TANX = TG(1)
      TANY = TG(2)
C
      DO M = 1,2
         DO I = 1,NPLA
C     
            NN = 0
            IF (M.EQ.2) NN = 1
            IF (MOD(I,2).EQ.NN) THEN
               SHIFT = +0.5
            ELSE
               SHIFT = -0.5
            ENDIF
C
            DISTY = -PIANO * (I - 1.)
            DISTX = -PIANO * (I - 1.)
C
            IF (M.EQ.1) THEN
               Y(I) = DISTX * TG(1) + CX
               BAR(M,I) = Y(I)
C
            ELSE
               YY(I) = DISTY * TG(2) + CY
               BAR(M,I) = YY(I)
C
            ENDIF
            CALL LASTRISCIA(BAR(M,I),IBAR(M,I))
        ENDDO
      ENDDO
C
      NCORE = 0.
      QCORE = 0.
      PLANEMAX = 1.01*(LOG(ABS(RIG)/0.0081)-1.)
      IPLANE = INT(ANINT(PLANEMAX)) + 5
      IF (IPLANE.GT.NPLA) IPLANE=NPLA
      DO J = 1,IPLANE
C
        NNX = IBAR(1,J)
        NNY = IBAR(2,J)
        IF (NNX.LT.9) NNX = 9
        IF (NNY.LT.9) NNY = 9
        IF (NNX.GT.88) NNX = 88
        IF (NNY.GT.88) NNY = 88
        INFX = NNX - 8
        INFY = NNY - 8
C
C 8 STRIPS ARE 2.88 cm , A MOLIERE RADIUS IS ABOUT 0.7 cm .
C
        ISUPX = NNX + 8
        ISUPY = NNY + 8
        RNSS = 0.
        QTOTT = 0.
        DO I = INFX,ISUPX
          IF (DEXY(1,J,I).GE.EMIN) THEN
            RNSS = RNSS + 1
            QTOTT = QTOTT + DEXY(1,J,I)
          ENDIF
        ENDDO
        DO I = INFY,ISUPY
          IF (DEXY(2,J,I).GE.EMIN) THEN
            RNSS = RNSS + 1
            QTOTT = QTOTT + DEXY(1,J,I)
          ENDIF
        ENDDO
        NCORE = RNSS * FLOAT(J) + NCORE
        QCORE = QTOTT * FLOAT(J) + QCORE
      ENDDO
C
      NINT = 0.
      CALL NOINT(NIN) ! if NINT=1 not interacting particle
      NINT = FLOAT(NIN)
C
C
      QCYL = 0.
      NCYL = 0.
C
C QCYL = DETECTED ENERGY AND NCYL = NUMBER OF HIT STRIPS IN A CYLINDER oF
C RADIUS 8.5 STRIPS WITH AXIS DEFINED BY THE DIRECTION OF THE INCOMING
C PARTICLE .
C
      DO J = 1,NPLA
C
        NNX = IBAR(1,J)
        NNY = IBAR(2,J)
        IF (NNX.LT.9) NNX = 9
        IF (NNY.LT.9) NNY = 9
        IF (NNX.GT.88) NNX = 88
        IF (NNY.GT.88) NNY = 88
        INFX = NNX - 8
        INFY = NNY - 8
C
C 8 STRIPS ARE 2.88 cm , A MOLIERE RADIUS IS ABOUT 0.7 cm .
C
        ISUPX = NNX + 8
        ISUPY = NNY + 8
        DO I = INFX,ISUPX
          IF (DEXY(1,J,I).LT.EMIN) GO TO 710
          NCYL = NCYL + 1
          QCYL = QCYL + DEXY(1,J,I)
  710   ENDDO
        DO I=INFY,ISUPY
          IF (DEXY(2,J,I).LT.EMIN) GO TO 810
          NCYL = NCYL + 1
          QCYL = QCYL + DEXY(1,J,I)
  810   ENDDO
      ENDDO
C
      QTRACK = 0.
      CALL LATERALE(QTRACK,RQT)
C
      QMAX = 0.
      DO M = 1,2
         DO I = 1,NPLA
            DO J = 1,NCHA
               IF (DEXY(M,I,J).GT.QMAX) QMAX = DEXY(M,I,J)
            ENDDO
         ENDDO
      ENDDO
C
      if (ival.eq.1) iev = iev + 1
      icount2 = icount2 + 1
      if (icount2.ge.nev) goto 50

c      print *,'qtot :',qtot,' nstrip :',nstrip,' ncore :',ncore,
c     &     ' qcore :',qcore

      call hfnt(1)
C
      goto 10

 50   continue

      print *,'Processed events :',icount2
      print *,'Events with no calculated baseline :',iev

***     /* Save histograms */

      CALL HROUT(0,icycle,' ')
      CALL HREND('Event')
      CLOSE (59)

      STOP
      END

C
C-------------------------------------------------------------------------
      SUBROUTINE CALCBASE(K,IA,IB,IPL,DX,BASE)

      IMPLICIT NONE

      INTEGER IA, IB, IPL, K
      REAL DX(11,96)
      REAL S1, S2
      REAL SIG
      REAL C1 
      REAL V(16)
      REAL RINF, RSUP
      REAL BASE, BAS
      REAL CHECK

      INTEGER J, I

      REAL PED(4,11,96),GOOD(4,11,96)

      COMMON / PEDE / PED, GOOD

      COMMON / CH / CHECK

C$$$      print *,'ia :',ia,' ib :',ib,' di :',di,' ds :',ds

      S1 = 0.
      S2 = 0.
      C1 = 0.
      DO J = IA,IB
         IF (GOOD(K,IPL,J).EQ.0) THEN 
            S1 = S1 + DX(IPL,J) 
            S2 = S2 + DX(IPL,J)*DX(IPL,J) 
            C1 = C1 + 1.
         ENDIF
      ENDDO
C
      BAS = 0.
      SIG = 0.
      IF (C1.GT.1) THEN
         BAS = S1 / C1
         SIG = SQRT(ABS((S2 - C1 * BAS*BAS) /(C1 - 1.)))
      ENDIF
C
      RINF = BAS - 3. * SIG
      RSUP = BAS + 3. * SIG
C
      S1 = 0.
      S2 = 0.
      C1 = 0.
      DO J = IA,IB
         IF (DX(IPL,J).GT.RINF.AND.DX(IPL,J).LT.RSUP
     &        .AND.GOOD(K,IPL,J).EQ.0) THEN
            S1 = S1 + DX(IPL,J) 
            S2 = S2 + DX(IPL,J)*DX(IPL,J) 
            C1 = C1 + 1.
         ENDIF
      ENDDO
C
      BAS = 0.
      SIG = 0.
      IF (C1.GT.1) THEN
         BAS = S1 / C1
         SIG = SQRT(ABS((S2 - C1 * BAS*BAS) /(C1 - 1.)))
      ENDIF
C
      IF (CHECK.EQ.1) PRINT *,'BAS :',BAS,' SIG :',SIG
      IF (SIG.LT.30) THEN
         BASE = BAS
         RETURN
      ENDIF
C
      RINF = 1.E6
      DO J = IA,IB
         IF (DX(IPL,J).GT.0.AND.DX(IPL,J).LT.RINF.AND.
     &        GOOD(K,IPL,J).EQ.0) THEN 
            RINF = DX(IPL,J)
            I = J
            IF (CHECK.EQ.1) PRINT *,'RINF :',RINF,' I :',I
         ENDIF
      ENDDO
C
      RINF = 1.E6
      DO J = IA,IB
         IF (DX(IPL,J).GT.0.AND.DX(IPL,J).LT.RINF.AND.
     &        GOOD(K,IPL,J).EQ.0.AND.J.NE.I) THEN
            RINF = DX(IPL,J)
            if (check.eq.1) print *,'rinf2 :',rinf
         ENDIF
      ENDDO
C
      S1 = 0.
      C1 = 0.
      DO J = IA,IB
         IF (DX(IPL,J).LT.(RINF+40.).AND.DX(IPL,J).GT.(RINF-40.)
     &        .AND.GOOD(K,IPL,J).EQ.0) THEN
            S1 = S1 + DX(IPL,J) 
            C1 = C1 + 1.
            if (check.eq.1) print *,'s1 :',s1,' c1 :',c1
         ENDIF   
      ENDDO

      IF (C1.GT.1) BASE = S1 / C1
      if (check.eq.1) print *,'BASE :',BASE
C      PRINT *,'BASE :',BASE,' C1 :',C1

      RETURN
      END
      
C
C---------------------------------------------------------------------
      SUBROUTINE LATERALE(RQT1,RQT2)
C---------------------------------------------------------------------
C RQT1 (IT WILL BE CALLED QTRACK IN THE N-TUPLE) IS THE SUM OF THE DETECTED
C ENERGY IN THE STRIP ALONG THE TRACK AND THE TWO CLOSEST STRIPS . FOR ALL THE
C LAYERS . RQT2 (IS NOT USED IN THE N-TUPLA) IS THE TOTAL ENERGY MINUS RQT1 .
C
      INCLUDE 'INTEST.TXT'
      REAL RQT1
      INTEGER A,B
      REAL BAR(2,NPLA)
      REAL Q(0:NPLA)
      INTEGER IBAR(2,NPLA)
      COMMON/ANGOLO/BAR,IBAR

      RQT2=0.

      INPIA = 1
C
      QQQ=0
      MAX=0
      Q(MAX)=0
C
      DO I = INPIA,NPLA
        A = IBAR(1,I)
        B = IBAR(2,I)
        IF (A.LE.2) A = 3
        IF (B.LE.2) B = 3
        IF (A.GE.(NCHA-1)) A = NCHA - 2
        IF (B.GE.(NCHA-1)) B = NCHA - 2

        DO J = A-1,A+1
          IF (DEXY(1,I,J).GE.EMIN) RQT1 = RQT1 + DEXY(1,I,J)
  600   ENDDO
C
        DO J = B-1,B+1
          IF (DEXY(2,I,J).GE.EMIN) RQT1 = RQT1 + DEXY(2,I,J)
        ENDDO
C
        DO J=1,A-2
          PXY = DEXY(1,I,J)
          IF (PXY.GE.EMIN) RQT2 = RQT2 + PXY
  650   ENDDO
C
        DO J=A+2,NCHA
          PXY = DEXY(1,I,J)
          IF (PXY.GE.EMIN) RQT2 = RQT2 + PXY
  700   ENDDO
C
        DO J=1,B-2
          PXY = DEXY(2,I,J)
          IF (PXY.GE.EMIN) RQT2 = RQT2 + PXY
  750   ENDDO
C
        DO J=B+2,NCHA
          PXY = DEXY(2,I,J)
          IF (PXY.GE.EMIN) RQT2 = RQT2 + PXY
  800   ENDDO
C
      ENDDO
C
C
  400 RETURN
      END

C------------------------------------------------
      SUBROUTINE PEDESTAL(ERROR)
C------------------------------------------------

      IMPLICIT NONE
C
      CHARACTER*6 file_name
C      
      integer icount,icount2
c
      integer*4 evnum, systime, len
      integer nev
      integer RUNERROR, ERROR
C
C Normal variables definition
C
      INTEGER FFD
C
      integer event
C
      INTEGER l, j, kk, m, mm, ll, ival, iev
      INTEGER i, istat, icycle, ierr, icheck, icheck2, icheck3
      INTEGER LVMIN

      integer nset,nword(10) 
C
      INTEGER*2 STRINGA(20000) 
C
      INTEGER ic, k, jj
      INTEGER status
      INTEGER inf, isup, idet

      INTEGER*2 length

      INTEGER lundata, iosop

      INTEGER*2 st, st2

      CHARACTER*2 nr
      CHARACTER*1 alfa

      REAL auto(7)

      REAL RP,RI

      REAL SUM(6), SUM2(6), COUNT(6), 
     &     BAS(4,11), DIS(6), DDIS(6)
      REAL SSUM, SSUM2, CCOUNT, MME, SSIG

      REAL VAL

      REAL DI, DS
      REAL ME(6), SIG(6)

      REAL PED(4,11,96), GOOD(4,11,96)
      COMMON / PEDE / PED, GOOD

C         
C Begin !
C
      RUNERROR=0                ! error variable

      PRINT *,'Pedestal file to read (only day and number included)? '
      READ(*,905)file_name

905   FORMAT(A6)

      lundata=44

      OPEN(UNIT=lundata,FILE='/wizard3/pamela/tb-0903/data/output_0309'
     +     //file_name//'.dat',STATUS='OLD',
     +  FORM='UNFORMATTED',ERR=50,IOSTAT=IOSOP)
C
      ffd=FNum(lundata)          ! take the file descriptor number
C

C Read event: use C readevent routine
C

 10   continue
C
      ival = 0
C
      call readev(evnum,systime,len,stringa,RUNERROR,ffd) 
      
      if (runerror.eq.-1.or.runerror.eq.1) goto 50

      ic = 0
 
      ic = ic + 1
      length = ic
      do kk = 1,4
         ic = length + 1
         k = ic
         write (*,41) stringa(k)
         st2 = IAND(stringa(k),'FF00'x)
         status = ISHFT(st2,-8)
         write (*,41) status
         status = ISHFT(st2,-8)
         IF (KK.eq.1.and.status.ne.170) then
            print *,'Problems with view:',KK
            ERROR = 1
            GOTO 50
         endif
         IF (KK.eq.2.and.status.ne.182) then
            print *,'Problems with view:',KK
            ERROR = 1
            GOTO 50
         endif
         IF (KK.eq.3.and.status.ne.177) then
            print *,'Problems with view:',KK
            ERROR = 1
            GOTO 50
         endif
         IF (KK.eq.4.and.status.ne.173) then
            print *,'Problems with view:',KK
            ERROR = 1
            GOTO 50
         endif
C
 41      format(1x,'status :',1x,Z4)
C
         ic = ic + 1
         length = length + (stringa(ic) + 2)
         print *,'Lenght1 :',length
         print *,'Lenght2 :',stringa(ic)

         if (stringa(ic).ne.4629) then
            print *,'problems with view',KK
            ERROR = 1
            GOTO 50
         endif
C
         do i = 1,11
            do j = 1,96
               ic = ic + 1
               ped(kk,i,j) = stringa(ic) 
               good(kk,i,j) = stringa(ic+1)
               if (kk.eq.1.and.i.eq.4.and.j.eq.77) 
     &              print *,'ped :',stringa(ic),' g/b :',
     %              stringa(ic+1)
               ic = ic + 1
            enddo
         enddo
C
      enddo
C

 50   continue


      CLOSE (44)

      RETURN
      END







