********************************************************
*
      SUBROUTINE CELLDIF(VPX,VPY,VPZ,ETAMAX,
     >                     SPX,SPY,SPZ,SE,FCALE,E_HK,
     >                     FEBP1,FEBP2)

********************************************************
* Routine to calculate various quantities 
* from calorimeter cells. The electron cells are removed.
* 
* NOTE: noise cut on calorimeter cells
*       needs to be discussed. Now is: 
* 	0.14 EM, 0.16  
*
* INPUT: vertex x,y,z = VP(1),VP(2),VP(3) 
*
*
*OUTPUT:
*
*       ZES name:      Local name:
*      ------------    ----------
*      ZES_Etamax_ce = ETAMAX
*      ZES_V_H_px_ce = SPX  	
*      ZES_V_H_py_ce = SPY  	
*      ZES_V_H_pz_ce = SPZ  
*      ZES_V_H_E_ce  = SE  	
*      ZES_E_FCAL    = FCALE
*      ZES_E_HK      = E_HK
*      ZES_Efcalbp_r1= FEBP1
*      ZES_Efcalbp_r2= FEBP2
*
******************************************************

      IMPLICIT NONE

#include "partap.inc"
#include "caltru.inc"
#include "sidat95.inc"
!  /*            ! for siNISTra  */

  
      INTEGER CAL, MODU,TOW,PMT,SEC,EH,I
      INTEGER NCELLS,ICAND,IEL,IERR
      REAL VPX,VPY,VPZ
      REAL FEBP1,FEBP2,FCALE
      REAL REMCE,BEMCE,FEMCE
      REAL ENER(3,32,23,2)
      REAL ETAC,E_HK,ETAMAX    
      REAL EMCCUT,HACCUT

      REAL VP(3),THETA_H,PHI_H,XP,YP,ZP
      REAL SPX,SPY,SPZ,SE
      LOGICAL FCALBP1,FCALBP2,FCALBP23,RCALBP1,RCALOUT
      LOGICAL CORNERS,FCALBP2N, FCALBP3
      LOGICAL RCALTOT,FCALTOT,BCALTOT
      LOGICAL RCALEMC,FCALEMC,BCALEMC
      INTEGER ELCELLS(100000)
      REAL       IMBCUT,ENECUT,IMBAL,ENERGY
      DATA       IMBCUT  / 0.90 /
      DATA       ENECUT  / 1.50 /


      DO CAL=1,3
	DO MODU=1,32
	  DO TOW=1,23
	    ENER(CAL,MODU,TOW,1)=0. ! EMC
	    ENER(CAL,MODU,TOW,2)=0. ! HAC
	  ENDDO
	ENDDO
      ENDDO
      CALL VZERO(ELCELLS,100000)
      VP(1)=VPX
      VP(2)=VPY
      VP(3)=VPZ


* Fill the array of the electron cells (Sinistra case)
* I assume Sinistra has already been called


      CALL FINDIS95( 1, 0.9, ICAND, IERR )   
      IF (ICAND.GT.0) THEN
        DO I = 1, POSDAT(3,ICAND)
          ELCELLS (POSDAT(I+3,ICAND)) = 1
        ENDDO
      ENDIF


      FEBP1 = 0.
      FEBP2 = 0.
      FCALE = 0.
      SPX   = 0.
      SPY   = 0.
      SPZ   = 0.
      SE    = 0.
      E_HK  = 0.
      ETAMAX = -100

* cut  on cell minimum energy
*          0.08 FOR EMC 
*          0.14 FOR HAC

      EMCCUT = 0.08
      HACCUT = 0.14

      DO I=1,COUTAB(CALTRU)
	CALL FETTAB(CALTRU,ID,I)


C Apply a standard cut on imbalance.
C === The logic is as follow: if it's a spark then the energy comes from only
C === one PMT and the imbalance is very high (90% and more)
C === For E<1.5 the imbalance can be very high.

        ENERGY = Caltru_E
        IF (CALTRU_imbal.EQ.0) THEN
          IMBAL = 0.
        ELSE
          IMBAL = Caltru_imbal/Caltru_E
        ENDIF
        IF (ABS(IMBAL) .GE. IMBCUT .AND. ENERGY.GT. ENECUT)  goto 123

C Remove the electron

        IF (ELCELLS(CALTRU_CELLNR).EQ.1) goto 123

	IF(CALTRU_E.GT.EMCCUT) THEN
	  PMT=CALTRU_CELLNR
	  CAL=0 ! must clear before call mvbits
	  CALL MVBITS(PMT,14,2,CAL,0)
	  CAL=CAL+1
	  MODU=0
	  CALL MVBITS(PMT,9,5,MODU,0)
	  MODU=MODU+1
	  TOW=0
	  CALL MVBITS(PMT,4,5,TOW,0)
	  TOW=TOW+1
	  SEC=0
	  CALL MVBITS(PMT,1,3,SEC,0)
	  IF(SEC.GT.4) EH=2  ! HAC
	  IF(SEC.LE.4) EH=1  ! EMC


          if(EH.EQ.2.AND.CALTRU_E.LT.HACCUT) GOTO 123


C ** first ring:module 11,12,13 and tower 11,12,13

	  FCALBP1=CAL.EQ.1.AND.ABS(MODU-12).LE.1.AND.
     &            ABS(TOW-12).LE.1

C ** second ring:module 10,11,12,13,14 and tower 10,11,12,13,14

	  FCALBP2=CAL.EQ.1.AND.ABS(MODU-12).LE.2.AND.
     &            ABS(TOW-12).LE.2

	  FCALTOT=CAL.EQ.1

	  IF(FCALBP1)FEBP1=FEBP1+CALTRU_E
	  IF(FCALBP2)FEBP2=FEBP2+CALTRU_E
	  IF(FCALTOT)FCALE=FCALE+CALTRU_E

* Calculate the angles to fill the components

          CALL ENGLE (CALTRU_CELLNR,VP,CALTRU_E,CALTRU_IMBAL,
     1                 THETA_H,PHI_H,XP,YP,ZP)

	   
           SPX    = SPX  + CALTRU_E*SIN(THETA_H)*COS(PHI_H)
           SPY    = SPY  + CALTRU_E*SIN(THETA_H)*SIN(PHI_H)
           SPZ    = SPZ  + CALTRU_E*COS(THETA_H)
           SE     = SE   + CALTRU_E

C ----     FILL ETAMAX AND ENERGY WITH ETA>3.1
C
           ETAC   = -LOG(TAN(THETA_H/2.))
           IF (ETAC.GE.3.1) THEN
             E_HK=E_HK+CALTRU_E
           ENDIF
           IF (ETAC.GE.ETAMAX) THEN	
             ETAMAX = ETAC
           ENDIF
	ENDIF

 123  ENDDO


      RETURN
      END