c-----------------------------------------------------------------
c userfunction for AYEBALL + PPAC + Neutrons + DSSD + Kitchen Sink
c
c Who's to blame: D.J. Blumenthal (standard parts)
c                 D. Seweryniak (tricky parts: neutrons & correlations)
c Last Updated: September 1995  
c-----------------------------------------------------------------

ccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccc
        logical function userin
ccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccc

c---
c  		User initialization routine
c		Called once at the beginning of aqcuisition and used
c		to initialize various things.
c---

	implicit none

	include 'aye_rdt_fix.inc'
	include 'aye_rdt_flt.inc'

c  		et 	- event type
c  		np	- number of parameters
c  		ev	- event vector
c
c..the pseudoevent structure is:
c
c	       	  pev(#)                 contents
c               -----------             -----------------------------
c		    1-20                Ege
c		    21-40               Tge
c		    41-60               doppler corrected energies 
c		    61-80               bit hit pattern germanium
c		    81-100              adc hit pattern germanium
c		    101,102,103	        1st, 2nd, 3rd corrected energy
c	      	    104		        number of eurogam Ge
c	      	    105		        number of tessa Ge
c	      	    106		        number of all Ge
c
c                   120                 event type bit register
c                   121                 event length
c
cxxxxx              122                 monitor DE
c                   123                 monitor E
c
c                   131                 R
c                   132                 L
c                   133                 U
c                   134                 D
c                   135                 X
c                   136                 RL               
c                   137                 Y
c                   138                 UD
c                   139                 deltaE
c                   140                 TDC TOF (PPAC - DSSD)
c                   141                 TAC TOF (PPAC - DSSD)
c                   142                 RF TOF 
c                   143                 rec-gamma TAC
c                   144                 PPAC pileup bit (1)
c                   145                 E * T**2
c                   146                 ppac DE in decay event
c
c                   200-215             n dets time of flight
c                   220-235             n dets total energy
c                   240-255             n dets slow component
c
c
c    302	adc hit bit register for implants
c    303        event type bit register for 2nd time around (fast decay)
c    304        adc hit bit register for decays
c
c    331     	emax			-recoil	-front
c    332	strip # with emax	-recoil	-front                  
c    333	total E of all strips	-recoil	-front
c    334        # of strips hit		-recoil	-front
c
c    335     	emax			-recoil	-back
c    336	strip # with emax	-recoil	-back                  
c    337	total E of all strips	-recoil	-back
c    338        # of strips hit		-recoil	-back
c
c    339        avg emax_rf & emax_rb
c    340	diff         "             (+2048)
c    341	avg etot_rf & etot_rb
c    342	diff         "             (+2048)
c
c    351     	emax			-decay	-front
c    352	strip # with emax	-decay	-front                  
c    353	total E of all strips	-decay	-front
c    354        # of strips hit		-decay	-front
c     
c    355     	emax			-decay	-back
c    356	strip # with emax	-decay	-back                  
c    357	total E of all strips	-decay	-back
c    358        # of strips hit		-decay	-back
c
c    359        avg emax_df & emax_db
c    360	diff         "             (+2048)
c    361	avg etot_df & etot_db
c    362	diff         "             (+2048)
c
c    365        (recoil front strip #) - (decay front strip #) + 100
c    366        (recoil  back strip #) - (decay  back strip #) + 100
c
c    311	time48 = time from 48 bit scaler in minutes (implant)
c    312          "                                          (decay)
c    313        (decay time) - (implant time) + 1000 when both in same event
c
c    320        focal plane Ge energy
c    321        focal plane Ge time    
c....correlation stuff
c
c    400        FMA x for a recoil followed by a decay
c    401        front energy for a decay correlated with a recoil
c    420        FMA x for a recoil followed by a decay 
c               time gated (t<t_lim_g1)
c    421        front energy for a decay correlated with a recoil
c               time gated (t<t_lim_g1)
c    402        back energy for a decay correlated with a recoil
c    403        time between recoil and a correlated decay in microseconds
c               divided by timecomp_c1
c    453        time between recoil and a correlated decay in microseconds
c               divided by timecomp_c2
c    463        time between recoil and a correlated decay in microseconds
c    404        front energy for a first decay 
c    405        front energy for a second decay 
c    424        front energy for a first decay 
c               time gated (t<t_lim_g1)
c    425        front energy for a second decay 
c               time gated (t<t_lim_g2)
c    406        time between first and second decay in microseconds
c               divided by timecomp_c1
c    456        time between first and second decay in microseconds
c               divided by timecomp_c2
c    466        time between first and second decay in microseconds
c    407        front energy for a second decay 
c    408        front energy for a third decay 
c    409        time between second and third decay in microseconds
c               divided by timecomp_c1
c    459        time between second and third decay in microseconds
c               divided by timecomp_c2
c    469        time between second and third decay in microseconds
c    410        correlation statistics:
c               1 - recoil
c               2 - decay
c               3 - recoil decay correltion
c               4 - decay decay correltion, first generation
c               5 - decay decay correltion, second generation
c               6 - decay decay correltion, third generation and higher

c-----------------------------------------------------------------------------
c..FIX:
c-----------------------------------------------------------------------------
c  1 - 48	OFF_RF(48)  -  offsets for gain matching FRONT recoil strips
c 49 - 96	OFF_RB(48)  -                            BACK
c 97 -144	OFF_DB(48)  -                            FRONT decay
c145 -192	OFF_DB(48)  -                            BACK
c   193 	ICAL        -  energy cal ON/OFF switch [1=on,0=off]
c   194         PULSER      -  flag = 0/1 if data is real/pulser
c   195         XET2GAIN    -  compression for E*T**2
c   196         corr_file
c   197         t_lim_g1
c   198         t_lim_g2
c   199         ermin       -  minimum channel # for recoil energy to be used
c-----------------------------------------------------------------------------
c.. FLT:
c-----------------------------------------------------------------------------                                                                
c  1 - 48	GAIN_RF(48) - gains for gain matching FRONT recoil strips
c 49 - 96	GAIN_RB(48) -                         BACK
c 97 -144	GAIN_DF(48) -                         FRONT decay
c145 -192	GAIN_DB(48) -                         BACK 
c   193 	XSCALE 	    - for PPAC X calculation  
c   194		XOFFSET     - for PPAC X calculation   
c   195		YSCALE      - for PPAC Y calculation
c   196		YOFFSET     - for PPAC Y calculation 
c   197 	TIMECOMP    - time compression so it fits in spectrum
c   198 	timecomp_c1
c   199 	timecomp_c2
c
c   231         BETA
c   232         GAIN 
c------------------------------------------------------------------------
c       declaring variables ...
c------------------------------------------------------------------------
	integer*2	et,np,type,size
	integer*2	ev(*),nPEV
	parameter	(nPEV=500)
	integer*2	pev(nPEV)

	integer*2	index,locev,locpev
        integer*4       i,j,k,m,n,ii
	integer*2	nrge1,nrge2,nrge
	integer*2	evtype,evbits(16),patge1(16),patge2(16)
        integer*2       npat(16)
        integer*4       nfold,nfoldn

        integer*2 fastdecay,evtype2
	real*8 time48,time48comp,time49,time49comp,timediff,timediff1
	real*4 stime0,stime1,stime2
	integer*4 cntbit
        integer*2 adchit(16)
	integer*2 stripmap_rf(16,3),stripmap_df(8,6)
	integer*2 stripmap_rb(16,3),stripmap_db(8,6)
	integer*2 strip,ch,nhit,adcid,teststrip
	integer*2 stripmax_rf,stripmax_rb,stripmax_df,stripmax_db
	integer*2 nhittot_rf,nhittot_rb,nhittot_df,nhittot_db
	real*4    e
	real*4 emax_rf,emax_rb, emax_df,emax_db
	real*4 etot_rf,etot_rb,etot_df,etot_db
	integer*2 bitmap_rf(3),bitmap_df(6), bitmap_rb(3),bitmap_db(6)   
c-------------------------------------------------------
        integer*2 rf_strip(48), rf_e(48), rb_strip(48), rb_e(48)
        integer*2 nhittot_rf_temp
        integer*2 nhittot_rb_temp,stripdiff
c-------------------------------------------------------      
                                                                   
	real*8		delran
	integer*4	iseed

	real*4		gamma
	real*4		ang1_5,dps1_5
	parameter	(ang1_5=158.)
	real*4		ang6_10,dps6_10
	parameter	(ang6_10=134.)
	real*4		ang11_15,dps11_15
	parameter	(ang11_15=101.0)
	real*4		ang16_20,dps16_20
	parameter	(ang16_20=79.0)
	real*4		degrad
	parameter	(degrad=57.29578)
	real*4		dps(20)

	real*4 R,L,U,D,X,Y,RL,UD
	integer*2 lowt,hight

        integer*2 data
        integer*4 ddata

        data lowt,hight /10,2040/
        data iseed /123457/

c-----------------------------------------------------------------        
	data stripmap_rf/ 
c..front mapping of adc chan (3 16ch adc's mapped to strips 1-48)
     1  1, 2, 3, 4, 5, 6, 7, 8, 9,10,11,12,13,14,15,16,
     2 17,18,19,20,21,22,23,24,25,26,27,28,29,30,31,32,
     3 33,34,35,36,37,38,39,40,41,42,43,44,45,46,47,48/     
	data stripmap_rb/ 
c..back mapping of adc chan (3 16ch adc's mapped to strips 1-48)
     1  1, 2, 3, 4, 5, 6, 7, 8, 9,10,11,12,13,14,15,16,
     2 17,18,19,20,21,22,23,24,25,26,27,28,29,30,31,32,
     3 33,34,35,36,37,38,39,40,41,42,43,44,45,46,47,48/     

	data stripmap_df/ 
c..front (5 8ch adc's mapped to strips 1-48)
     1  1, 2, 3, 4, 5, 6, 7, 8, 
     2  9,10,11,12,13,14,15,16,
     3 17,18,19,20,21,22,23,24,
     4 25,26,27,28,29,30,31,32,
     5 33,34,35,36,37,38,39,40,
     6 41,42,43,44,45,46,47,48/          

	data stripmap_db/ 
c..back (5 8ch adc's mapped to strips 1-48)
     1  1, 2, 3, 4, 5, 6, 7, 8, 
     2  9,10,11,12,13,14,15,16,
     3 17,18,19,20,21,22,23,24,
     4 25,26,27,28,29,30,31,32,
     5 33,34,35,36,37,38,39,40,
     6 41,42,43,44,45,46,47,48/          

c..BiRa input register has a bit set for each of the 6 recoil adc's
c  12,11,10 (9,8,7) are for the front (back).
	data bitmap_rf/12,11,10/
	data bitmap_rb/9,8,7/
c..bit exists for decay adc's in pairs only so 2 adc's should be read
c  out for each decay adc bit which is set.
	data bitmap_df/12,12,11,11,10,10/
	data bitmap_db/9,9,8,8,7,7/

        integer*4 spn, mass, d1, d12, d123

        logical user, user1, user2, user3, user4, userex, fixflt

        INTEGER*4  INZVECTOR(*),W2VECTOR(*),W1VECTOR(*),CONVECTOR(*)

        CHARACTER*10 CON_NAME, H_NAME
        LOGICAL H_OK, CON_OK 
        LOGICAL CONNUMBER
        LOGICAL H1NUMBER, H2NUMBER
        INTEGER*4 Z1_CON_N, Z2_CON_N, Z3_CON_N, Z4_CON_N, Z5_CON_N, Z6_CON_N
        INTEGER*4 M1_CON_N, M2_CON_N, M3_CON_N, M4_CON_N, M5_CON_N
        INTEGER*4 N_CON_N(0:15)
        INTEGER*4 NE_H_N, NSE_H_N, NRF_H_N
        INTEGER*4 NHIT_H_N, NHITN_H_N, NHITG_H_N
        INTEGER*4 NHITRF_H_N, NHITE_H_N, NHITSE_H_N
        INTEGER*4 GGMN_H_N, GMN_H_N
        INTEGER*4 D1_N(0:4),D12_N(0:4),D123_N(0:4)
        INTEGER*4 rgtac_n
        INTEGER*4 GGDSSDN_H_N, GDSSDN_H_N
        INTEGER*4 hitge,hitge2,ege_spec(20),tge_spec(20)
c----------------
c correlations
c----------------
        CHARACTER*6 REC_CON_NAME
        LOGICAL REC_CON_N_OK
        INTEGER*4 REC_CON_N
        CHARACTER*6 dec_CON_NAME
        LOGICAL dec_CON_N_OK
        INTEGER*4 dec_CON_N

        integer*2 dssd(0:17,48,48)
        real*8 dssd_time(48,48)

        character*9 date_str,time_str
        integer ioerr

        real dec_en_comp_corr /200.0/
        real rec_en_comp_corr /2000.0/
        real x_comp_corr /4.0/
        real time_comp_corr /1000.0/
        real tof_comp_corr /16.0/
                        
        character*40 corr_file_name /'corr_file.dat'/
        character*5 e_unit /' '/
        character*4 t_unit /'(ms)'/
        character*4 x_unit /'(/4)'/

c----------------------------------------------------------------------
c        getting access to conditions ...
c----------------------------------------------------------------------
        CON_NAME='M1'
        CON_OK=CONNUMBER(CON_NAME,M1_CON_N)
        CON_NAME='M2'
        CON_OK=CONNUMBER(CON_NAME,M2_CON_N)
        CON_NAME='M3'
        CON_OK=CONNUMBER(CON_NAME,M3_CON_N)
        CON_NAME='M4'
        CON_OK=CONNUMBER(CON_NAME,M4_CON_N)
        CON_NAME='M5'
        CON_OK=CONNUMBER(CON_NAME,M5_CON_N)

        CON_NAME='N0'
        CON_OK=CONNUMBER(CON_NAME,N_CON_N(0))
        CON_NAME='N1'
        CON_OK=CONNUMBER(CON_NAME,N_CON_N(1))
        CON_NAME='N2'
        CON_OK=CONNUMBER(CON_NAME,N_CON_N(2))
        CON_NAME='N3'
        CON_OK=CONNUMBER(CON_NAME,N_CON_N(3))
        CON_NAME='N4'
        CON_OK=CONNUMBER(CON_NAME,N_CON_N(4))
        CON_NAME='N5'
        CON_OK=CONNUMBER(CON_NAME,N_CON_N(5))
        CON_NAME='N6'
        CON_OK=CONNUMBER(CON_NAME,N_CON_N(6))
        CON_NAME='N7'
        CON_OK=CONNUMBER(CON_NAME,N_CON_N(7))
        CON_NAME='N8'
        CON_OK=CONNUMBER(CON_NAME,N_CON_N(8))
        CON_NAME='N9'
        CON_OK=CONNUMBER(CON_NAME,N_CON_N(9))
        CON_NAME='N10'
        CON_OK=CONNUMBER(CON_NAME,N_CON_N(10))
        CON_NAME='N11'
        CON_OK=CONNUMBER(CON_NAME,N_CON_N(11))
        CON_NAME='N12'
        CON_OK=CONNUMBER(CON_NAME,N_CON_N(12))
        CON_NAME='N13'
        CON_OK=CONNUMBER(CON_NAME,N_CON_N(13))
        CON_NAME='N14'
        CON_OK=CONNUMBER(CON_NAME,N_CON_N(14))
        CON_NAME='N15'
        CON_OK=CONNUMBER(CON_NAME,N_CON_N(15))

        REC_CON_NAME='REC'
        REC_CON_N_OK=CONNUMBER(REC_CON_NAME,REC_CON_N)

        DEC_CON_NAME='DEC'
        DEC_CON_N_OK=CONNUMBER(DEC_CON_NAME,DEC_CON_N)         

        CON_NAME='D1_1'
        CON_OK=CONNUMBER(CON_NAME,D1_N(0))

        CON_NAME='D1_2'
        CON_OK=CONNUMBER(CON_NAME,D1_N(1))

        CON_NAME='D1_3'
        CON_OK=CONNUMBER(CON_NAME,D1_N(2))

        CON_NAME='D1_4'
        CON_OK=CONNUMBER(CON_NAME,D1_N(3))

        CON_NAME='D1_5'
        CON_OK=CONNUMBER(CON_NAME,D1_N(4))

        CON_NAME='D12_1'
        CON_OK=CONNUMBER(CON_NAME,D12_N(0))

        CON_NAME='D12_2'
        CON_OK=CONNUMBER(CON_NAME,D12_N(1))

        CON_NAME='D12_3'
        CON_OK=CONNUMBER(CON_NAME,D12_N(2))

        CON_NAME='D12_4'
        CON_OK=CONNUMBER(CON_NAME,D12_N(3))

        CON_NAME='D12_5'
        CON_OK=CONNUMBER(CON_NAME,D12_N(4))

        CON_NAME='D123_1'
        CON_OK=CONNUMBER(CON_NAME,D123_N(0))

        CON_NAME='D123_2'
        CON_OK=CONNUMBER(CON_NAME,D123_N(1))

        CON_NAME='D123_3'
        CON_OK=CONNUMBER(CON_NAME,D123_N(2))

        CON_NAME='D123_4'
        CON_OK=CONNUMBER(CON_NAME,D123_N(3))

        CON_NAME='D123_5'
        CON_OK=CONNUMBER(CON_NAME,D123_N(4))

        CON_NAME='RGTAC'
        CON_OK=CONNUMBER(CON_NAME,RGTAC_N)

c        
c        getting access to histograms ...
c

        H_NAME='NE'
        H_OK=H2NUMBER(H_NAME,NE_H_N)
        H_NAME='NSE'
        H_OK=H2NUMBER(H_NAME,NSE_H_N)
        H_NAME='NRF'
        H_OK=H2NUMBER(H_NAME,NRF_H_N)

        H_NAME='NHIT'
        H_OK=H1NUMBER(H_NAME,NHIT_H_N)

        H_NAME='NHITRF'
        H_OK=H1NUMBER(H_NAME,NHITRF_H_N)
        H_NAME='NHITE'
        H_OK=H1NUMBER(H_NAME,NHITE_H_N)
        H_NAME='NHITSE'
        H_OK=H1NUMBER(H_NAME,NHITSE_H_N)

        H_NAME='NHITN'
        H_OK=H1NUMBER(H_NAME,NHITN_H_N)
        H_NAME='NHITG'
        H_OK=H1NUMBER(H_NAME,NHITG_H_N) 

        H_NAME='GGMN'
        H_OK=H2NUMBER(H_NAME,GGMN_H_N)

        H_NAME='GMN'
        H_OK=H2NUMBER(H_NAME,GMN_H_N)

        H_NAME='GGDSSDN'
        H_OK=H2NUMBER(H_NAME,GGDSSDN_H_N)

        H_NAME='GDSSDN'
        H_OK=H2NUMBER(H_NAME,GDSSDN_H_N)


c..germanium spectra
        H_NAME='HITGE'
        H_OK=H1NUMBER(H_NAME,HITGE)
        H_NAME='HITGE2'
        H_OK=H1NUMBER(H_NAME,HITGE2)

        do i=1,20
           if (i.lt.10) then
             write(h_name,'(''E'',i1) ') i
           else
             write(h_name,'(''E'',i2) ') i
           endif
           H_OK=H1NUMBER(H_NAME,ege_spec(i))
           if (i.lt.10) then
             write(h_name,'(''T'',i1) ') i
           else
             write(h_name,'(''T'',i2) ') i
           endif
           H_OK=H1NUMBER(H_NAME,tge_spec(i))
        enddo
c------------------------------------------------------------------

	call defPEV (pev(1),pev(nPEV))
	call deffix(fixarr,fixzzz)
	call defflt(fltarr,fltzzz)
c---
c		Define the random number seeds.
c---
	iseed = 123457

        userin=.true.

c        
c       opening a correlation file ...
c

        if (corr_file.eq.1) then
          write (*,*) 'OPENING A CORRELATION FILE ...'
          open(unit=10,file=CORR_FILE_name,status='unknown',
     @         access='append',shared,iostat=ioerr,err=100)
          write(10,*) '********************************************************************************'
          call date(date_str)
          call time(time_str) 
          write(10,*)  'OPEN FILE'
          write(10,*) date_str//' - '//time_str
          write(10,*) '********************************************************************************'
        endif

        return

100     write (*,*) 'CANNOT OPEN CORRELATION FILE !'

ccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccc
	entry  user (et,np,ev,INZVECTOR,W2VECTOR,W1VECTOR,CONVECTOR)
ccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccc

c-------------------------------------------------------------
c  START PROCESSING THE EVENT HERE!!!!!!!! 
c------------------------------------------------------------- 
        user=.false.
        call zeroI2(npev,pev)
	call zeroI2(16,adchit)
	call zeroI2(16,evbits)
	evtype=0
	evtype2=0
        fastdecay=0
        time48=0
        time49=0
        timediff=0         
	nrge1 = 0
	nrge2 = 0
        nrge  = 0      
        locev=1                       !pointer for raw event vector
        pev(121)=np                   !event length in words

c-------------------------------------------------------------
c  get the event type bit register 
c------------------------------------------------------------- 
        evtype=ev(locev)
        locev=locev+1
        call bit16(evtype,evbits)
        pev(120)=evtype

        if (evbits(8).eq.1 .and. evbits(9).eq.1) then
cxxxx           type*,'#########IMP+DEC------>kill decay part'
           evbits(9)=0
        endif

c-------------------------------------------------------------
c  Monitor Event  
c-------------------------------------------------------------
      if ( evbits(1) .eq. 1) then

        pev(123) = and (ev(locev),'0FFF'x) !emon
        locev=locev+1
        userin=.true.
        return                               !ignore anything else
      endif  

c-------------------------------------------------------------
c	IMPLANTATION
c-------------------------------------------------------------
	   if ( evbits(8).eq.1 ) then
c-----
c time 
c-----
   	      stime0=ev(locev)
	      stime1=ev(locev+1)
	      stime2=ev(locev+2)
	      locev=locev+3
       	      if (stime0.lt.0) then
	         stime0 = stime0 + 65536
	      endif
              if (stime1.lt.0) then
	         stime1 = stime1 + 65536
	       endif
              if (stime2.lt.0) then
                 stime2 = stime2 + 65536
	      endif
  	      time48 = stime2 * (2 ** 16)
	      time48 = ( time48 + stime1 ) * (2 ** 16)
	      time48 = time48 + stime0
	      time48comp = time48 / timecomp
	      if (time48comp.lt.32768) then
	         pev(311) = int(time48comp)
	      endif
c-------------------------
c ADC hit pattern register
c-------------------------
	   call bit16(ev(locev),adchit)
	   pev(302)=ev(locev)
           locev=locev+1

           locpev=331
c---------------------------
c recoil adc readout - FRONT
c---------------------------
c..i=1-3 ---> adc #
c..j=1-? ---> # of chan read in a given adc
	         e=0.
	         emax_rf=0.
	         stripmax_rf=0
	         etot_rf=0.
	         nhittot_rf=0
                 if (pulser.eq.1) then
                    if (teststrip.ge.48) teststrip=0
                    teststrip=teststrip+1
                 endif
	         do i=1,3
c..adc 1-3 are front, 4-6 are back so must use convoluted index for adchit
	            if (adchit(bitmap_rf(i)).eq.1) then
		       nhit=cntbit(ev(locev))
	               locev=locev + 1
		       do j=1,nhit
    		          ch = ishft(ev(locev),-12) + 1
		          strip = stripmap_rf(ch,i)
    		          e = iand(ev(locev),'0FFF'X)
		          if (ical.eq.1) then
		             e = (e+ran(iseed)-0.5)*gain_rf(strip)
     1                                                 +off_rf(strip)
		          endif
                          if (e.gt.ermin) then
                             nhittot_rf=nhittot_rf+1
                             rf_strip(nhittot_rf)=strip
                             rf_e(nhittot_rf)=e
                          endif
                          if (pulser.eq.1) then
                             if (strip.eq.teststrip) then
		                stripmax_rf=strip
	         	        emax_rf=e
                             endif
                          endif
		          locev=locev + 1
                       enddo   
                    endif
                enddo    

                if (pulser.eq.0) then
c..compute corrected fold & strip energies (put neighboring hits together)
                nhittot_rf_temp = nhittot_rf
                if (nhittot_rf.gt.1) then
                  do i=1,nhittot_rf-1
                    stripdiff=rf_strip(i+1)-rf_strip(i)
                    if (stripdiff.eq.1 .or. stripdiff.eq.-1) then
                        nhittot_rf=nhittot_rf-1
                        if (rf_e(i+1).gt.rf_e(i)) then
                          rf_e(i+1)=rf_e(i)+rf_e(i+1)                          
                          rf_e(i)=0
                        else
                          rf_e(i)=rf_e(i)+rf_e(i+1)
                          rf_e(i+1)=0
                        endif
                    endif 
                  enddo   
                endif
c..compute emax,etot,strip # with emax
                do i=1,nhittot_rf_temp
                  if (rf_e(i).gt.emax_rf) then
                    emax_rf=rf_e(i)
                    stripmax_rf=rf_strip(i)                    
                  endif
                  etot_rf=etot_rf+rf_e(i)
                enddo   
              
                endif

c..max energy
	         pev(locpev)=nint(emax_rf)
c..strip # containing max energy
	         pev(locpev+1)=stripmax_rf
c..total energy in all strips
                 if (etot_rf.gt.32768) etot_rf=32768
                 pev(locpev+2)=nint(etot_rf)
c..# of strips which fired
	         pev(locpev+3)=nhittot_rf

	         locpev=locpev+4
c---------------------------
c recoil adc readout - BACK
c---------------------------
	         e=0.
	         emax_rb=0.
	         stripmax_rb=0
	         etot_rb=0.
	         nhittot_rb=0
	         do i=1,3
	            if (adchit(bitmap_rb(i)).eq.1) then
		       nhit=cntbit(ev(locev))
	               locev=locev + 1
		       do j=1,nhit
    		          ch = ishft(ev(locev),-12) + 1
		          strip = stripmap_rb(ch,i)
    		          e = iand(ev(locev),'0FFF'X)
		          if (ical.eq.1) then
		             e = (e+ran(iseed)-0.5)*gain_rb(strip)
     1                                                 +off_rb(strip)
		          endif
                          if (e.gt.ermin) then
                             nhittot_rb=nhittot_rb+1
                             rb_strip(nhittot_rb)=strip
                             rb_e(nhittot_rb)=e
                          endif
                          if (pulser.eq.1) then
                             if (strip.eq.teststrip) then
		                stripmax_rb=strip
	         	        emax_rb=e
                             endif
	                  endif
		          locev=locev + 1
                       enddo    
                    endif
                enddo

                 if (pulser.eq.0) then
c..compute corrected fold & strip energies (put neighboring hits together)
                 nhittot_rb_temp = nhittot_rb
                 if (nhittot_rb.gt.1) then
                   do i=1,nhittot_rb-1
                    stripdiff=rb_strip(i+1)-rb_strip(i)
                    if (stripdiff.eq.1 .or. stripdiff.eq.-1) then
                         nhittot_rb=nhittot_rb-1
                         if (rb_e(i+1).gt.rb_e(i)) then
                           rb_e(i+1)=rb_e(i)+rb_e(i+1)                          
                           rb_e(i)=0
                         else
                           rb_e(i)=rb_e(i)+rb_e(i+1)
                           rb_e(i+1)=0
                         endif
                     endif 
                   enddo   
                 endif
c..compute emax,etot,strip # with emax
                 do i=1,nhittot_rb_temp
                   if (rb_e(i).gt.emax_rb) then
                     emax_rb=rb_e(i)
                     stripmax_rb=rb_strip(i)                    
                   endif
                   etot_rb=etot_rb+rb_e(i)
                 enddo    
 
                 endif

c..max energy
	         pev(locpev)=nint(emax_rb)
c..strip # containing max energy
	         pev(locpev+1)=stripmax_rb
c..total energy in all strips
                 if (etot_rb.gt.32768) etot_rb=32768
	         pev(locpev+2)=nint(etot_rb)
c..# of strips which fired
	         pev(locpev+3)=nhittot_rb

	         locpev=locpev+4
c------------------------------------------------------
c recoil F/B averages, differences of emax's and totals
c------------------------------------------------------
c..avg emax_r
	      pev(339) = nint( (emax_rf + emax_rb) * 0.5 )
c..f/b diff in emax_r
	      pev(340) = emax_rf - emax_rb + 2048
c..avg etot_r
	      pev(341) = nint( (etot_rf + etot_rb) * 0.5 )
c..f/b diff in etot_r
	      pev(342) = etot_rf - etot_rb + 2048

	      locpev=343

           endif  
c----------^^^^^--------------------------
c end of implant loop (JUST THE DSSD PART)
c-----------------------------------------

c-----------------------
c  Do all the PPAC stuff
c-----------------------
      if ( iand(evtype,'00AA'x) .ne. 0) then     !ppac,rg,rgg,implant

c..first generate the (X,Y) position & the LR and UD sums        
              R = FLOAT(iand(EV(locev),'0FFF'X))
              L = FLOAT(iand(EV(locev+1),'0FFF'X))
              U = FLOAT(iand(EV(locev+2),'0FFF'X))
              D = FLOAT(iand(EV(locev+3),'0FFF'X))
c..store R,L,U,D readouts (raw data)
              PEV(131)=R
              PEV(132)=L
              PEV(133)=U
              PEV(134)=D
              IF(R.GT.lowt .AND .R.LT. hight
     &                     .AND.
     &           L.GT.lowt .AND. L.LT. hight) THEN 
                 X  = ( (R-L) + XOFFSET ) * XSCALE   
                 RL = (R+L)/2
                 PEV(135)=NINT(X)
                 PEV(136)=NINT(RL)
              ENDIF        
              IF(U.GT.lowt .AND. U.LT. hight
     &                     .AND.
     &           D.GT.lowt .AND.D.LT. hight) THEN
                 Y = ( (D-U) + YOFFSET ) * YSCALE   
                 UD = (U+D)/2
                 PEV(137)=NINT(Y)
                 PEV(138)=NINT(UD)
             ENDIF

    	     locev = locev + 4
                
c..PPAC DELTA ENERGY
          pev(139)= and (ev(locev),'0FFF'x)
      	  locev = locev + 1
c..TOF tdc (ppac-dssd time of flight)   
             pev(140) = and (ev(locev),'0FFF'x)
	     locev = locev + 1
c..TOF TAC (ppac-dssd time of flight for long times)   
             pev(141) = and (ev(locev),'0FFF'x)
	     locev = locev + 1
c..TOF (ppac-rf)
	     pev(142) = and (ev(locev),'0FFF'x)
	     locev = locev + 1
c..TOF (rec-gamma TAC)
	     pev(143) = and (ev(locev),'0FFF'x)
	     locev = locev + 1
c...pileup bit for PPAC
           pev(144)= and (ev(locev),'0FFF'x)
           LOCEV=LOCEV+1   
c..E * T**2
c            if (pev(162).gt.1 .and. pev(140).gt.1) then
c               pev(145) = nint(float(pev(162)) 
c     1                 * float(pev(140))**2 / float(xet2gain))
c            endif    


      endif    !end of ppac block

c-----------------------------------------------------------
c  Do all the germanium and neutron stuff here if they fired
c-----------------------------------------------------------
      if ( iand(evtype,'00BC'x) .ne. 0) then    !g,rg,gg,rg,implant

c..Get the hit pattern for the germaniums.
	call bit16(ev(locev),patge1)
	call bit16(ev(locev+1),patge2)
        locev=locev+2

 	do i=1,20
	   j = i       !loc of raw ge e in pev
	   k = 20 + i  !loc of ge time in pev
	   l = 40 + i  !loc of corrected ge e in pev
	   m = 60 + i  !loc of bit hit patt for ge #i in pev    
	   n = 80 + i  !loc of adc hit patt for ge #i in pev

           index=0
           if (i.le.10) then
              if (patge1(i).eq.1) index=1                !EG are 1-10
           else                                         
              if (patge2(i-10).eq.1) index=2             !TESSA are 11-20
           endif                           
	   if (index.gt.0) then                          !check ge bits
              pev(j) = ev(locev)                         !raw ge energy
              ddata=pev(j)
              call h1inc(ege_spec(i),ddata)
              pev(k) = iand(ev(locev+1),'0fff'x)         !ge time
              ddata=nint(float(pev(k))/2.0)
              call h1inc(tge_spec(i),ddata)
              locev = locev + 2
              call h1inc(HITGE,i)                        !bit hit pattern
              if (pev(j).gt.1 .and. pev(j).lt.8190) then !valid ge energy
                 call h1inc(HITGE2,i)                    !adc hit pattern
	         nrge = nrge + 1                         !total ge mult
                 if (index.eq.1) then
   	            nrge1 = nrge1 + 1                    !EG mult
                 else
                    nrge2 = nrge2 + 1                    !TESSA mult
                 endif
	         delran = ran(iseed) - 0.5
                 pev(l) = nint( (dps(i)*pev(j) + delran)/2.0 )
                 if ( i.ne.9 .and. i.ne.10) then         !exclude LEPS from 1st 3 gams
                    if (nrge.le.3) pev(100+nrge)=pev(l)  !1st 3 energies
                 endif
              endif 
	   endif
        enddo    

	pev(104) = nrge1
	pev(105) = nrge2
        pev(106) = nrge	

c-------------------------------------
c   take care of the neutron detectors
c-------------------------------------
        call bit16(ev(locev),npat)        
            
        locev=locev+1

        locpev=200

        nfold=0
 	do i=0,15
          if ((npat(i+1)).eq.1) then

c           neutron detector hit pattern       

            II=I
            call h1inc(NHIT_H_N,ii)

c           neutron detector RF       

            data=iand(ev(locev),'0FFF'X)
            pev(locpev+i)=data
            ddata=data
            call h2inc(NRF_H_N,ddata,ii)   

c           neutron detector RF hit pattern      
            if ((data.gt.0).and.(data.lt.4095)) call h1inc(NHITRF_H_N,ii)

c           neutron detector energy                     

            data=iand(ev(locev+1),'0FFF'X)
            pev(locpev+20+i)=data
            ddata=data
            call h2inc(NE_H_N,ddata,ii)   

c           neutron detector energy hit pattern      
            if ((data.gt.0).and.(data.lt.2047)) call h1inc(NHITE_H_N,ii)

c           neutron detector slow component energy                     

            data=iand(ev(locev+2),'0FFF'X)
            pev(locpev+40+i)=data
            ddata=data
            call h2inc(NSE_H_N,ddata,ii)   

c           neutron detector energy hit pattern      
            if ((data.gt.0).and.(data.lt.2047)) call h1inc(NHITSE_H_N,ii)

            locev=locev+3
            nfold=nfold+1
          endif
        enddo

        pev(260)=nfold

      endif         !end of germanium/neutron block

c--------------------------------------------------------------------------
c if the ppac and gamma stuff were at the end of an implantation event, now
c determine if normal recoil event or one with following fast decay event
c--------------------------------------------------------------------------
      if (evbits(8).eq.1) then  
         evtype2 = ev(locev)                     !read evtype again
         pev(303)=evtype2
         locev=locev+1
         if (iand(evtype2,256) .ne. 0) then   !if fast decay we wrote out a 256
c            type*,'fastdecay, evtype2= ',evtype2
            fastdecay = 1
         endif	                              
      endif

c.......................................................................	
c	DECAY event
c.......................................................................
      if ( evbits(9).eq.1 .or. fastdecay.eq.1) then
         if (fastdecay.ne.1) then
            locev=2
         endif
c.....
c time
c.....
 	      stime0=ev(locev)
	      stime1=ev(locev+1)
	      stime2=ev(locev+2)
              locev=locev+3     
	      if (stime0.lt.0) then
	         stime0 = stime0 + 65536
	      endif
                if (stime1.lt.0) then
	           stime1 = stime1 + 65536
	         endif
                if (stime2.lt.0) then
                   stime2 = stime2 + 65536
	      endif
   	      time49 = stime2 * (2 ** 16)
 	      time49 = ( time49 + stime1 ) * (2 ** 16)
 	      time49 = time49 + stime0
 	      time49comp = time49 / timecomp
 	      if (time49comp.lt.32768) then
 	         pev(312) = int(time49comp)
 	      endif

c------------------
c PPAC DELTA ENERGY
c------------------
              pev(146)= iand (ev(locev),'0FFF'x)
              locev=locev+1

c-------------------------
c ADC hit pattern register
c-------------------------
	      call bit16(ev(locev),adchit)
	      pev(304)=ev(locev)
              locev=locev+1

   	      locpev=351

c------------------------
c decay adc readout - FRONT
c------------------------
	         e=0.
	         emax_df=0.
	         stripmax_df=0
	         etot_df=0.
	         nhittot_df=0
                 if (pulser.eq.1) then
                    if (teststrip.ge.48) teststrip=0
                    teststrip=teststrip+1
                 endif
	         do i=1,6
	            if (adchit(bitmap_df(i)).eq.1) then
c..extract adc id # we defined
	             adcid=iand(ev(locev),'00FF'x)
c..check which adc in pair fired. If it is the I th one, do stuff.
	             if (adcid.eq.i) then
c..get nhit from header word
		       nhit=and(rshift(ev(locev),8),'000F'x)
	               nhittot_df=nhittot_df +nhit
	               locev=locev + 2
		       do j=1,nhit
c..NOTE: mask for silena adc's has 7, not 15 since only 3 bits are chan id
    		          ch=iand(ishft(ev(locev),-12),'0007'x)+1
		          strip = stripmap_df(ch,i)
    		          e = iand(ev(locev),'0FFF'X)
		          if (ical.eq.1) then
		             e = (e+ran(iseed)-0.5)*gain_df(strip)
     1                                                 +off_df(strip)
		          endif
                          if (pulser.eq.0) then
		             if (e.gt.emax_df) then
                                stripmax_df=strip
                                emax_df=e
                             endif
                          else
                             if (strip.eq.teststrip) then
		                stripmax_df=strip
	         	        emax_df=e
                             endif
	                  endif
ccccccccccccccc		          etot_df=etot_df + e
		          locev=locev + 1
                       enddo    
	             endif
                    endif
                 enddo      
c..max energy
	         pev(locpev)=nint(emax_df)
c..strip # containing max energy
	         pev(locpev+1)=stripmax_df
c..total energy in all strips
                 if (etot_df.gt.32768) etot_df=32768
	         pev(locpev+2)=nint(etot_df)
c..# of strips which fired
	         pev(locpev+3)=nhittot_df

	         locpev=locpev+4
c-------------------------
c decay adc readout - BACK
c-------------------------
	         e=0.
	         emax_db=0.
	         stripmax_db=0
	         etot_db=0.
	         nhittot_db=0
	         do i=1,6
	            if (adchit(bitmap_db(i)).eq.1) then
	             adcid=and(ev(locev),'00FF'x)
	             if ( adcid .eq. (i+6) ) then
		       nhit=iand(rshift(ev(locev),8),'000F'x)
	               nhittot_db=nhittot_db +nhit
	               locev=locev + 2
		       do j=1,nhit
    		          ch=iand(ishft(ev(locev),-12),'0007'x)+1
		          strip = stripmap_db(ch,i)
    		          e = iand(ev(locev),'0FFF'X)
		          if (ical.eq.1) then
		             e = (e+ran(iseed)-0.5)*gain_db(strip)
     1                                                 + off_db(strip)
		          endif
                          if (pulser.eq.0) then
		             if (e.gt.emax_db) then
                                stripmax_db=strip
                                emax_db=e
                             endif
                          else
                             if (strip.eq.teststrip) then
		                stripmax_db=strip
	         	        emax_db=e
                             endif
	                  endif
cccccccccccccc		          etot_db=etot_db + e
		          locev=locev + 1
                       enddo   
	             endif
                    endif
                 enddo   
c---max energy
	         pev(locpev)=nint(emax_db)
c---strip # containing max energy
	         pev(locpev+1)=stripmax_db
c---total energy in all strips
                 if (etot_db.gt.32768) etot_db=32768
	         pev(locpev+2)=nint(etot_db)
c---number of strips which fired
	         pev(locpev+3)=nhittot_db

	         locpev=locpev+4
c------------------------------------------------------
c decay F/B averages, differences of emax's and totals
c------------------------------------------------------
c..avg emax_d
	      pev(359) = nint( (emax_df + emax_db) * 0.5 )
c..f/b diff in emax_d
	      pev(360) = emax_df - emax_db + 2048
c..avg etot_d
	      pev(361) = nint( (etot_df + etot_db) * 0.5)
c..f/b diff in etot_d
	      pev(362) = etot_df - etot_db + 2048
c-------------------------------------------------------------
c time difference between implant and decay when in same event
c-------------------------------------------------------------
              if (fastdecay.eq.1) then
                 pev(365)=pev(332)-pev(352) + 100
                 pev(366)=pev(336)-pev(356) + 100
                 if (time48.gt.0 .and. time49.gt.0 
     1               .and. time49.gt.time48) then
                    timediff=time49-time48 + 1000.
                 endif
                 if (timediff.lt.32768) then              
                    pev(312)=int(timediff)
                 endif
              endif

c---------------
c focal plane ge
c---------------
             pev(320)=ev(locev)        !E
             pev(321)=ev(locev+1)      !t
             locev=locev+2

      endif
c-----^^^^^--------
c end of decay loop
c------------------

      user=.true.

      return

c=======================================================================
C     MASS and NEUTRON FOLd
ccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccc
      entry user1(et,np,ev,INZVECTOR,W2VECTOR,W1VECTOR,CONVECTOR)
ccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccc
c     
c     getting mass ...
c

      mass=0
      if (convector(2*m1_con_n-1).eq.1) then
        mass=1
      elseif (convector(2*m2_con_n-1).eq.1) then
        mass=2
      elseif (convector(2*m3_con_n-1).eq.1) then
        mass=3         
      elseif (convector(2*m4_con_n-1).eq.1) then
        mass=4
      elseif (convector(2*m5_con_n-1).eq.1) then
        mass=5
      endif
c     
c     getting neutron fold ...
c
   
      nfoldn=0
      do i=0,15
        if (convector(2*N_CON_n(i)-1).eq.1) then
          nfoldn=nfoldn+1
          ii=i
          call h1inc(nhitn_h_n,ii)
        endif
      enddo

      pev(261)=nfoldn

c
c     incrementing gamma spectra gated with Z, M/Q and neutron fold ...
c
      spn=mass*24+nfoldn

      if (nrge.eq.1) then
        ii=pev(101)
        call h2inc(gmn_h_n,ii,spn)
      else
        do i=1,nrge
          ii=pev(100+i)
          call h2inc(ggmn_h_n,ii,spn)
        enddo    
      endif

      user1=.true.

      return

c=======================================================================
C       CORRELATION PART FOR A RECOIL EVENT
c=======================================================================

c
c       dssd matrix:
c
c       0 - pixel flag
c       1 - PPAC X
c       2 - emax_rf
c       3 - ppac-dssd time of flight
c       4 - emax_df #1
c       5 - emax_db #1
c       6 - emax_df #2
c       7 - emax_db #2
c       8 - emax_df #3
c       9 - emax_db #3

c      10 - t1/time_comp_g1
c      11 - t2/time_comp_g2
c      12 - t3/time_comp_g3

c      13 - number of ge detectors
c      14 - number of ge detectors
c      15 - number of ge detectors
c      16 - number of ge detectors

c       pseudo events:

c       400 - PPAC X
c       401 - emax_df #1
c       402 - emax_db #1
c       404 - emax_df #1 (first generation)  
c       405 - emax_df #2 
c       407 - emax_df #2 (second generation) 
c       408 - emax_df #3

c       410 - correlation status

c       430 - PPAC X (defined for first generation)
c       431 - t1/time_comp_g1 (defined for first generation) 
c       432 - PPAC X (defined for second generation)
c       433 - t1/time_comp_g1 (defined for second generation) 
c       434 - emax_df #1 (defined for second generation) 
c       435 - t2/time_comp_g2 (defined for second generation)

c       403 - t1/timecomp_c1
c       453 - t1/time_comp_c2
c       463 - t1
c       473 - t1/time_comp_g1
c       406 - t2/timecomp_c1
c       456 - t2/time_comp_c2
c       466 - t2
c       476 - t2/time_comp_g2
c       409 - t3/timecomp_c1
c       459 - t3/time_comp_c2
c       469 - t3
c       479 - t3/time_comp_g2

c=======================================================================

	entry user2(type,size,ev,inzvector,w2vector,w1vector,convector)        

c=======================================================================

c        type *,'user2'        

        user2=.false.

        if (evbits(8).eq.1) then

C        SELECT RECOILS BY CHECKING 'REC' GATE

          if (convector(2*rec_con_n-1).eq.1) then
            user2=.true.
            dssd(0,stripmax_rf,stripmax_rb)=1
            dssd(1,stripmax_rf,stripmax_rb)=pev(135)
            dssd(2,stripmax_rf,stripmax_rb)=emax_rf
            dssd(3,stripmax_rf,stripmax_rb)=pev(140)
            dssd_time(stripmax_rf,stripmax_rb)=time48
            pev(410)=1                             
            ! now gamma rays !
            dssd(13,stripmax_rf,stripmax_rb)=nrge 
            do i=1,3
              dssd(13+i,stripmax_rf,stripmax_rb)=pev(100+i)
            enddo
            dssd(17,stripmax_rf,stripmax_rb)=pev(261) 
          endif
        endif        

        return

c=======================================================================
C       CORRELATION PART FOR A DECAY EVENT
c=======================================================================

	entry user3(type,size,ev,inzvector,w2vector,w1vector,convector)        

c=======================================================================

c        type *,'user3'

        user3=.false.

        if ((evbits(9).eq.1).or.(fastdecay.eq.1)) then

C         SELECT DECAYS BY CHECKING 'DEC' GATE

          if (convector(2*dec_con_n-1).eq.1) then

          user3=.true.
          pev(410)=2

          if (dssd(0,stripmax_df,stripmax_db).eq.1) then

C           RECOIL-DECAY CORRELATION

            dssd(0,stripmax_df,stripmax_db)=2
            pev(410)=3
            dssd(4,stripmax_df,stripmax_db)=emax_df 
            dssd(5,stripmax_df,stripmax_db)=emax_db 
            pev(400)=dssd(1,stripmax_df,stripmax_db)
            pev(401)=emax_df 
            pev(402)=emax_db 
            timediff1=(time49-dssd_time(stripmax_df,stripmax_db))
     @               /timecomp_c1            
            if ((timediff1.lt.256).and.(timediff1.gt.0)) then
              pev(403)=iidint(timediff1)
            else
              pev(403)=255
            endif
            timediff1=(time49-dssd_time(stripmax_df,stripmax_db))
     @              /timecomp_c2
            if ((timediff1.lt.256).and.(timediff1.gt.0)) then
              pev(453)=iidint(timediff1)
            else
              pev(453)=255
            endif               
            timediff1=time49-dssd_time(stripmax_df,stripmax_db)
            if ((timediff1.lt.256).and.(timediff1.gt.0)) then
              pev(463)=iidint(timediff1)
            else
              pev(463)=255
            endif 

            timediff1=(time49-dssd_time(stripmax_df,stripmax_db))
     @      /timecomp_g1
            if ((timediff1.lt.256).and.(timediff1.gt.0)) then
              pev(473)=iidint(timediff1)
            else
              pev(473)=255
            endif 

            dssd(10,stripmax_df,stripmax_db)=pev(473)

c            pev(464)=iidint(timediff1)              

c            if ((timediff1.lt.t_lim_g1).and.(timediff1.gt.0)) then
c              pev(420)=dssd(1,stripmax_df,stripmax_db)
c              pev(421)=emax_df 
c           endif

c           if (timediff1.le.500) then
c             pev(497)=dssd(2,stripmax_df,stripmax_db)
c           endif

           timediff1=time49-dssd_time(stripmax_df,stripmax_db)

           dssd_time(stripmax_df,stripmax_db)=time49
            if (corr_file.eq.1) then
              write(10,1000) stripmax_df,stripmax_db, 
     @        e_unit,
     @        emax_df/dec_en_comp_corr,emax_db/dec_en_comp_corr,
     @        x_unit,
     @        nint(dssd(1,stripmax_df,stripmax_db)/x_comp_corr),
     @        dssd(2,stripmax_df,stripmax_db)/dec_en_comp_corr,
     @        nint(dssd(3,stripmax_df,stripmax_db)/tof_comp_corr),
     @        t_unit,
     @        timediff1/time_comp_corr

1000          format('RD0>',' F-B',2i3,' EF-EB',a5,2f6.3,
     @        ' X',a4,i4,' ER',f6.2,'TOF',i4,' T',a4,f10.3)
            endif
          elseif (dssd(0,stripmax_df,stripmax_db).eq.2) then

C             DECAY-DECAY CORRELATION FIRST GENERATION

              dssd(0,stripmax_df,stripmax_db)=3
              pev(410)=4
              dssd(6,stripmax_df,stripmax_db)=emax_df ! or etot_df
              dssd(7,stripmax_df,stripmax_db)=emax_db ! or etot_df
              pev(405)=emax_df
              pev(404)=dssd(4,stripmax_df,stripmax_db)
              timediff1=(time49-dssd_time(stripmax_df,stripmax_db))
     @                 /timecomp_c1
              if ((timediff1.lt.256).and.(timediff1.gt.0)) then
                pev(406)=iidint(timediff1)
              else
                pev(406)=255
              endif
              timediff1=(time49-dssd_time(stripmax_df,stripmax_db))
     @                 /timecomp_c2
              if ((timediff1.lt.256).and.(timediff1.gt.0)) then
                pev(456)=iidint(timediff1)
              else
                pev(456)=255
              endif               
              timediff1=time49-dssd_time(stripmax_df,stripmax_db)
              if ((timediff1.lt.256).and.(timediff1.gt.0)) then
                pev(466)=iidint(timediff1)
              else
                pev(466)=255
              endif

              timediff1=(time49-dssd_time(stripmax_df,stripmax_db))
     @        /timecomp_g2
              if ((timediff1.lt.256).and.(timediff1.gt.0)) then
                pev(476)=iidint(timediff1)
              else
                pev(476)=255
              endif 

              dssd(11,stripmax_df,stripmax_db)=pev(476)

c              if ((timediff1.lt.t_lim_g2).and.(timediff1.gt.0)) then
c                pev(425)=emax_df
c                pev(424)=dssd(4,stripmax_df,stripmax_db)
c              endif

              pev(430)=dssd(1,stripmax_df,stripmax_db)               
              pev(431)=dssd(10,stripmax_df,stripmax_db)               

              timediff1=time49-dssd_time(stripmax_df,stripmax_db)
               
              dssd_time(stripmax_df,stripmax_db)=time49

              if (corr_file.eq.1) then
                write(10,1020) stripmax_df,stripmax_db, 
     @          e_unit,
     @          emax_df/dec_en_comp_corr,emax_db/dec_en_comp_corr,
     @          e_unit,
     @          float(dssd(4,stripmax_df,stripmax_db))/dec_en_comp_corr, 
     @          x_unit,
     @          nint(dssd(1,stripmax_df,stripmax_db)/x_comp_corr),
     @          t_unit,
     @          timediff1/time_comp_corr



1020            format('DD1>',' F-B',2i3,' EF2-EB2',a5,2f6.3,
     @          ' EF1',a5,f6.3,' X',a4,i4,
     @          ' T',a4,f10.3)

              endif
              
              elseif (dssd(0,stripmax_df,stripmax_db).eq.3) then

C             DECAY-DECAY CORRELATION SECOND GENERATION

              dssd(0,stripmax_df,stripmax_db)=4
              pev(410)=5
              dssd(8,stripmax_df,stripmax_db)=emax_df ! or etot_df
              dssd(9,stripmax_df,stripmax_db)=emax_db ! or etot_db

              pev(407)=dssd(6,stripmax_df,stripmax_db)
              pev(408)=emax_df ! or etot_df
              timediff1=(time49-dssd_time(stripmax_df,stripmax_db))
     @                 /timecomp_c1
              if ((timediff1.lt.256).and.(timediff1.gt.0)) then
                pev(409)=iidint(timediff1)
              else
                pev(409)=255
              endif
              timediff1=(time49-dssd_time(stripmax_df,stripmax_db))
     @                 /timecomp_c2
              if ((timediff1.lt.256).and.(timediff1.gt.0)) then
                pev(459)=iidint(timediff1)
              else
                pev(459)=255
              endif               
              timediff1=time49-dssd_time(stripmax_df,stripmax_db)
              if ((timediff1.lt.256).and.(timediff1.gt.0)) then
                pev(469)=iidint(timediff1)
              else
                pev(469)=255
              endif               

              timediff1=(time49-dssd_time(stripmax_df,stripmax_db))
     @        /timecomp_g3
              if ((timediff1.lt.256).and.(timediff1.gt.0)) then
                pev(479)=iidint(timediff1)
              else
                pev(479)=255
              endif 

              dssd(12,stripmax_df,stripmax_db)=pev(479)

              pev(432)=dssd(1,stripmax_df,stripmax_db)               
              pev(433)=dssd(10,stripmax_df,stripmax_db)               
              pev(434)=dssd(4,stripmax_df,stripmax_db)               
              pev(435)=dssd(11,stripmax_df,stripmax_db)               

              timediff1=time49-dssd_time(stripmax_df,stripmax_db)

              dssd_time(stripmax_df,stripmax_db)=time49

              if (corr_file.eq.1) then          

                write(10,1040) stripmax_df,stripmax_db, 
     @          e_unit,
     @          emax_df/dec_en_comp_corr,emax_db/dec_en_comp_corr,
     @          e_unit,
     @          float(dssd(6,stripmax_df,stripmax_db))/dec_en_comp_corr, 
     @          e_unit,
     @          float(dssd(4,stripmax_df,stripmax_db))/dec_en_comp_corr, 
     @          x_unit,
     @          nint(dssd(1,stripmax_df,stripmax_db)/x_comp_corr),
     @          t_unit,
     @          timediff1/time_comp_corr


1040            format('DD2>',' F-B',2i3,' EF3-EB3',a5,2f6.3,
     @          ' EF2',a5,f6.3,' EF1',a5,f6.3,
     @          ' X',a4,i4,' T',a4,f10.3)

              endif

              elseif (dssd(0,stripmax_df,stripmax_db).gt.3) then

C             DECAY-DECAY CORRELATION HIGHER THAN SECOND GENERATION

              pev(410)=6
              if (corr_file.eq.1) then          
                write(10,*) 'DD3> not included !'
              endif 
          endif                 
        endif
        endif

        return

c=======================================================================
c       DSSD gates       
c=======================================================================
	entry user4(type,size,ev,inzvector,w2vector,w1vector,convector)        
c=======================================================================

c     
c     checking DSSD conditions ...
c

      user4=.false.

      if (convector(2*rgtac_n-1).eq.1) then

      d1=0
      d12=0
      d123=0

      if (convector(2*d1_n(0)-1).eq.1) then
        d1=1
      elseif (convector(2*d1_n(1)-1).eq.1) then
        d1=2
      elseif (convector(2*d1_n(2)-1).eq.1) then
        d1=3         
      elseif (convector(2*d1_n(3)-1).eq.1) then
        d1=4
      elseif (convector(2*d1_n(4)-1).eq.1) then
        d1=5
      endif

      if (convector(2*d12_n(0)-1).eq.1) then
        d12=1
      elseif (convector(2*d12_n(1)-1).eq.1) then
        d12=2
      elseif (convector(2*d12_n(2)-1).eq.1) then
        d12=3         
      elseif (convector(2*d12_n(3)-1).eq.1) then
        d12=4
      elseif (convector(2*d12_n(4)-1).eq.1) then
        d12=5
      endif

      if (convector(2*d123_n(0)-1).eq.1) then
        d123=1
      elseif (convector(2*d123_n(1)-1).eq.1) then
        d123=2
      elseif (convector(2*d123_n(2)-1).eq.1) then
        d123=3         
      elseif (convector(2*d123_n(3)-1).eq.1) then
        d123=4
      elseif (convector(2*d123_n(4)-1).eq.1) then
        d123=5
      endif

c
c     incrementing gamma spectra gated with DSSD and neutron fold ...
c

      nfoldn=dssd(17,stripmax_df,stripmax_db)
      nrge=dssd(13,stripmax_df,stripmax_db)

      if (nfoldn.gt.3) nfoldn=3 
      spn=d1*4+nfoldn
      if (nrge.eq.1) then
        ii=dssd(14,stripmax_df,stripmax_db)
        call h2inc(gdssdn_h_n,ii,spn)
      else
        do i=1,nrge
          ii=dssd(13+i,stripmax_df,stripmax_db)
          call h2inc(ggdssdn_h_n,ii,spn)
        enddo    
      endif

      spn=d12*4+nfoldn+24
      if (nrge.eq.1) then
        ii=dssd(14,stripmax_df,stripmax_db)
        call h2inc(gdssdn_h_n,ii,spn)
      else
        do i=1,nrge
          ii=dssd(13+i,stripmax_df,stripmax_db)
          call h2inc(ggdssdn_h_n,ii,spn)
        enddo    
      endif

      spn=d123*4+nfoldn+48
      if (nrge.eq.1) then
        ii=dssd(14,stripmax_df,stripmax_db)
        call h2inc(gdssdn_h_n,ii,spn)
      else
        do i=1,nrge
          ii=dssd(13+i,stripmax_df,stripmax_db)
          call h2inc(ggdssdn_h_n,ii,spn)
        enddo    
      endif

      user4=.true.

      endif

      return

c=======================================================================
ccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccc
	entry userex
ccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccc
        if (corr_file.eq.1) then
          write(*,*) 'CLOSING A CORRELATION FILE ...'
          write(10,*) '********************************************************************************'
          call date(date_str)
          call time(time_str) 
          write(10,*)  'CLOSE FILE'
          write(10,*) date_str//' - '//time_str
          write(10,*) '********************************************************************************'
          close(unit=10)
        endif

        userex=.true.

	return

ccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccc
	entry fixflt
ccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccc
c---
c		Entry FixFlt is called each time the fix/flt parameters
c		are being changed. In this way the doppler shifts always
c		have the appropriate values
c---

	gamma = 1./sqrt(1. - beta*beta)

	dps1_5 = (1. - beta*cos(ang1_5/degrad))*gamma
	dps6_10 = (1. - beta*cos(ang6_10/degrad))*gamma
	dps11_15 = (1. - beta*cos(ang11_15/degrad))*gamma
	dps16_20 = (1. - beta*cos(ang16_20/degrad))*gamma
	
	dps(1)  = dps1_5    
	dps(2)  = dps1_5      
	dps(3)  = dps1_5 
	dps(4)  = dps1_5 
	dps(5)  = dps1_5 
	dps(6)  = dps6_10    
	dps(7)  = dps6_10
	dps(8)  = dps6_10    
	dps(9)  = dps6_10
	dps(10) = dps6_10
	dps(11) = dps11_15
	dps(12) = dps11_15   
	dps(13) = dps11_15
	dps(14) = dps11_15
	dps(15) = dps11_15 
	dps(16) = dps16_20
	dps(17) = dps16_20
	dps(18) = dps16_20 
	dps(19) = dps16_20 
	dps(20) = dps16_20

        fixflt=.true.

	return
	end