ToFLevel2/src/toftrk.for

*****************************************************************************
      INTEGER FUNCTION TOFTRK()

C****************************************************************************
C     31-08-06  WM
C     Changed to use DOTRACK2
C     Beta calculation: now the flightpath (instead of cos(theta)) is used
C     Beta calculation: all 4 TDV measurements must be < 4095 (in the old
C     routine it was (t1+t2)<8000
C
C     08-12-06 WM:
C     adc_c-bug :  The raw ADC value was multiplied with cos(theta)
C     and AFTER that there was an if statement "if tof32(right,i,iadc) < 4095"
C
C  jan-07 GF: ADC/TDCflags(4,12) inserted to flag artificial ADC/TDC
C             values
C  jan-07 WM: artificial ADC values created using attenuation calibration
C  jan-07 WM: artificial TDC values created using xy_coor calibration
C  jan-07 WM: modified xtofpos flag "101". xtofpos must be inside physical
C             dimension of the paddle +/- 10 cm
C  jan-07 WM: if xtofpos=101 then this paddle is not used for beta
C             calculation
C  jan-07 WM: in the xtofpos calculation a check for TDC.ne.4095 was
C             inserted. In the old code one would still calculate a
C             xtofpos-value even if the TDC information was missing
C  jan-07 WM: flag for PMTs #10 and #35 added, TDC=819 due to bit-shift
C  jan-05 WM: bug fixed: calculation of zenith angles using DOTRACK2
C             was incorrect
C  jan-07 WM: bug fixed: in some cases tdc_tw was calculated due to a
C             leftover "xhelp" value
C  apr-07 WM: attenuation fit curve is now a double exponential fit
C             conversion from raw ADC to pC using calibration function
C             new variables xtr_tof(6) and ytr_tof(6) give track position
C             at ToF layers
C  aug-07 WM: artificial ADC creation revised: Now an ADC value is created
C             only if there is a TDC value (before ADC was created in ANY
C             case)
C  jan-08 WM: Major Update: Time Walk correction introduced
C             Additionally we use the information from the "check_charge" 
C             function to fill artificial ADC values and make small corrections
C             to the k1-parameter (for Z>2)
C  feb-08 WM: Calculation of beta(13) changed: First a mean beta is calculated,
C             then in a second step we check the residuals of the single
C             measurements, reject if > 10 sigma, calculate chi2 and "quality"
C             beta is taken as good if chi2<20 and quality>10
C             The function "newbeta" is located in "tofl2com.for"
C  mar-08 WM: Call to "newbeta" changed, now a flag tells the function if the
C             call comes from "tofl2com" or form "toftrack"
C  mar-08 WM: Bug found in dEdx if check_charge>1
C  apr-08 WM: Bug found in S22 artificial ADC, mismatch found between the track
C             length  from DOTRACK2 and "GetLength" method for 4 combinations
C  oct-08 WM: New method to create artificial ADC values. Do NOT take the position
C             from the tracking, but the position from timing. This method gives a
C             better time resolution
C  nov-09 WM: the dEdx part ("adc_c") moved to the new dEdx routine from Napoli
C
C****************************************************************************
      IMPLICIT NONE
C
      include  'input_tof.txt'
      include  'output_tof.txt'
      include  'tofcomm.txt'
C

c     =======================================
c     variables for tracking routine
c     =======================================
      integer NPOINT_MAX
      parameter(NPOINT_MAX=100)

c     define TOF Z-coordinates
      integer NPTOF
      parameter (NPTOF=6)
      DOUBLE PRECISION ZTOF(NPTOF)
      DATA ZTOF/53.74,53.04,23.94,23.44,-23.49,-24.34/ !Sergio 9.05.2006

      integer itof,pmt_id

      DOUBLE PRECISION al_p(5),
     &     xout(NPOINT_MAX),yout(NPOINT_MAX),zin(NPTOF),
     &     THXOUT(NPOINT_MAX),THYOUT(NPOINT_MAX),TLOUT(NPOINT_MAX)


      INTEGER IFAIL
c      REAL dx,dy,dr
      REAL ds
      REAL t1,t2,t3,t4
      REAL yhelp,xhelp,xhelp1,xhelp2
      REAL c1,c2
C     REAL sw,sxw,w_i
      REAL dist,dl,F
      INTEGER ievent
C      INTEGER icount
C      REAL beta_mean
      REAL btemp(12)
      REAL hepratio

      INTEGER j,hitvec(6)

      real atten,pc_adc,check_charge,newbeta


      REAL theta,phi
C--   DATA ZTOF/53.74,53.04,23.94,23.44,-23.49,-24.34/ !Sergio 9.05.2006
      REAL tofarm12
      PARAMETER (tofarm12 = 29.70) ! from 53.39 to 23.69
      REAL tofarm23
      PARAMETER (tofarm23 = 47.61) ! from 23.69 to -23.92
      REAL tofarm13
      PARAMETER (tofarm13 = 77.31) ! from 53.39 to -23.92


C--- new
      REAL xtofpre(3),ytofpre(3)

      REAL y_coor_lin11c(8,2),x_coor_lin12c(6,2)
      REAL x_coor_lin21c(2,2),y_coor_lin22c(2,2)
      REAL y_coor_lin31c(3,2),x_coor_lin32c(3,2)

      DATA y_coor_lin11c(1,1),y_coor_lin11c(1,2) /-20.66,-2.497/
      DATA y_coor_lin11c(2,1),y_coor_lin11c(2,2) /-9.10, -2.52/
      DATA y_coor_lin11c(3,1),y_coor_lin11c(3,2) /-24.07,-2.12/
      DATA y_coor_lin11c(4,1),y_coor_lin11c(4,2) /-13.40,-2.47/
      DATA y_coor_lin11c(5,1),y_coor_lin11c(5,2) /-31.07,-2.32/
      DATA y_coor_lin11c(6,1),y_coor_lin11c(6,2) /-21.69,-2.63/
      DATA y_coor_lin11c(7,1),y_coor_lin11c(7,2) /-12.37,-2.65/
      DATA y_coor_lin11c(8,1),y_coor_lin11c(8,2) /-10.81,-3.15/

      DATA x_coor_lin12c(1,1),x_coor_lin12c(1,2) /12.96, -2.65/
      DATA x_coor_lin12c(2,1),x_coor_lin12c(2,2) /17.12,-2.44/
      DATA x_coor_lin12c(3,1),x_coor_lin12c(3,2) /7.26, -1.98/
      DATA x_coor_lin12c(4,1),x_coor_lin12c(4,2) /-22.52,-2.27/
      DATA x_coor_lin12c(5,1),x_coor_lin12c(5,2) /-18.54,-2.28/
      DATA x_coor_lin12c(6,1),x_coor_lin12c(6,2) /-7.67,-2.15/

      DATA x_coor_lin21c(1,1),x_coor_lin21c(1,2) /22.56,-1.56/
      DATA x_coor_lin21c(2,1),x_coor_lin21c(2,2) /13.94,-1.56/

      DATA y_coor_lin22c(1,1),y_coor_lin22c(1,2) /-24.24,-2.23/
      DATA y_coor_lin22c(2,1),y_coor_lin22c(2,2) /-45.99,-1.68/

      DATA y_coor_lin31c(1,1),y_coor_lin31c(1,2) /-22.99,-3.54/
      DATA y_coor_lin31c(2,1),y_coor_lin31c(2,2) /-42.28,-4.10/
      DATA y_coor_lin31c(3,1),y_coor_lin31c(3,2) /-41.29,-3.69/

      DATA x_coor_lin32c(1,1),x_coor_lin32c(1,2) /0.961, -3.22/
      DATA x_coor_lin32c(2,1),x_coor_lin32c(2,2) /4.98,-3.48/
      DATA x_coor_lin32c(3,1),x_coor_lin32c(3,2) /-22.08,-3.37/

C---

      INTEGER ihelp
      REAL xkorr,xpos

      INTEGER IZ
      REAL k1corrA1,k1corrB1,k1corrC1

      REAL yl,yh,xl,xh
C
      REAL hmemor(9000000)
      INTEGER Iquest(100)
C
      DATA ievent / 0 /

      INTEGER ifst
      DATA ifst /0/

      COMMON / pawcd / hmemor
      save / pawcd /
C
      Common / QUESTd / Iquest
      save / questd /
C
C     Begin !
C
      TOFTRK = 0

*******************************************************************

       if (ifst.eq.0) then
       ifst=1

C  ratio helium to proton ca. 4
      hepratio = 4.

      offset = 1
      slope = 2
      left = 1
      right = 2
      none_ev = 0
      none_find = 0
      tdc_ev = 1
      adc_ev = 1
      itdc = 1
      iadc = 2


      k1corrA1 = 0.
      k1corrB1 = -5.0
      k1corrC1=  8.0

      ENDIF   ! ifst

*******************************************************************

      ievent = ievent +1

      do i=1,13
         beta_a(i) = 100.
      enddo

      do i=1,4
         do j=1,12
            adc_c(i,j) = 1000.
         enddo
      enddo

      do i=1,12
         do j=1,4
            tofmask(j,i) = 0
         enddo
      enddo

      do i=1,4
         do j=1,12
            adcflag(i,j) = 0
         enddo
      enddo

      do i=1,4
         do j=1,12
            tdcflag(i,j) = 0
         enddo
      enddo

      pmt_id=0

      do j=1,6
      THXOUT(j) = 0.
      THYOUT(j) = 0.
      enddo

      do j=1,6
      xtr_tof(j) = 100.
      ytr_tof(j) = 100.
      enddo


C----------------------------------------------------------------------
C-------------------------get ToF data --------------------------------
C     we cannot use the tofxx(x,x,x)  data  from tofl2com since it is
C     manipulated (Time-walk, artificila ADc and TDC values using ToF
C     standalone information
C----------------------------------------------------------------------

c     put the adc and tdc values from ntuple into tofxx(i,j,k) variables

      do j=1,8
         tof11(1,j,2) = pc_adc(adc(ch11a(j),hb11a(j)))
         tof11(2,j,2) = pc_adc(adc(ch11b(j),hb11b(j)))
         tof11(1,j,1) = (tdc(ch11a(j),hb11a(j)))
         tof11(2,j,1) = (tdc(ch11b(j),hb11b(j)))
      enddo


      do j=1,6
         tof12(1,j,2) = pc_adc(adc(ch12a(j),hb12a(j)))
         tof12(2,j,2) = pc_adc(adc(ch12b(j),hb12b(j)))
         tof12(1,j,1) = (tdc(ch12a(j),hb12a(j)))
         tof12(2,j,1) = (tdc(ch12b(j),hb12b(j)))
      enddo

      do j=1,2
         tof21(1,j,2) = pc_adc(adc(ch21a(j),hb21a(j)))
         tof21(2,j,2) = pc_adc(adc(ch21b(j),hb21b(j)))
         tof21(1,j,1) = (tdc(ch21a(j),hb21a(j)))
         tof21(2,j,1) = (tdc(ch21b(j),hb21b(j)))
      enddo
 
      do j=1,2
         tof22(1,j,2) = pc_adc(adc(ch22a(j),hb22a(j)))
         tof22(2,j,2) = pc_adc(adc(ch22b(j),hb22b(j)))
         tof22(1,j,1) = (tdc(ch22a(j),hb22a(j)))
         tof22(2,j,1) = (tdc(ch22b(j),hb22b(j)))
      enddo

      do j=1,3
         tof31(1,j,2) = pc_adc(adc(ch31a(j),hb31a(j)))
         tof31(2,j,2) = pc_adc(adc(ch31b(j),hb31b(j)))
         tof31(1,j,1) = (tdc(ch31a(j),hb31a(j)))
         tof31(2,j,1) = (tdc(ch31b(j),hb31b(j)))
      enddo

      do j=1,3
         tof32(1,j,2) = pc_adc(adc(ch32a(j),hb32a(j)))
         tof32(2,j,2) = pc_adc(adc(ch32b(j),hb32b(j)))
         tof32(1,j,1) = (tdc(ch32a(j),hb32a(j)))
         tof32(2,j,1) = (tdc(ch32b(j),hb32b(j)))
      enddo

C----------------------------------------------------------------------

      DO i = 1,8
         if (abs(tof11(1,i,itdc)).gt.10000.) tof11(1,i,itdc)= 10000.
         if (abs(tof11(2,i,itdc)).gt.10000.) tof11(2,i,itdc)= 10000.
         if (abs(tof11(1,i,iadc)).gt.10000.) tof11(1,i,iadc)= 10000.
         if (abs(tof11(2,i,iadc)).gt.10000.) tof11(2,i,iadc)= 10000.
      ENDDO

      DO i = 1,6
         if (abs(tof12(1,i,itdc)).gt.10000.) tof12(1,i,itdc)= 10000.
         if (abs(tof12(2,i,itdc)).gt.10000.) tof12(2,i,itdc)= 10000.
         if (abs(tof12(1,i,iadc)).gt.10000.) tof12(1,i,iadc)= 10000.
         if (abs(tof12(2,i,iadc)).gt.10000.) tof12(2,i,iadc)= 10000.
      ENDDO


      DO i = 1,2
         if (abs(tof21(1,i,itdc)).gt.10000.) tof21(1,i,itdc)= 10000.
         if (abs(tof21(2,i,itdc)).gt.10000.) tof21(2,i,itdc)= 10000.
         if (abs(tof21(1,i,iadc)).gt.10000.) tof21(1,i,iadc)= 10000.
         if (abs(tof21(2,i,iadc)).gt.10000.) tof21(2,i,iadc)= 10000.
      ENDDO

      DO i = 1,2
         if (abs(tof22(1,i,itdc)).gt.10000.) tof22(1,i,itdc)= 10000.
         if (abs(tof22(2,i,itdc)).gt.10000.) tof22(2,i,itdc)= 10000.
         if (abs(tof22(1,i,iadc)).gt.10000.) tof22(1,i,iadc)= 10000.
         if (abs(tof22(2,i,iadc)).gt.10000.) tof22(2,i,iadc)= 10000.
      ENDDO

      DO i = 1,3
         if (abs(tof31(1,i,itdc)).gt.10000.) tof31(1,i,itdc)= 10000.
         if (abs(tof31(2,i,itdc)).gt.10000.) tof31(2,i,itdc)= 10000.
         if (abs(tof31(1,i,iadc)).gt.10000.) tof31(1,i,iadc)= 10000.
         if (abs(tof31(2,i,iadc)).gt.10000.) tof31(2,i,iadc)= 10000.
      ENDDO

      DO i = 1,3
         if (abs(tof32(1,i,itdc)).gt.10000.) tof32(1,i,itdc)= 10000.
         if (abs(tof32(2,i,itdc)).gt.10000.) tof32(2,i,itdc)= 10000.
         if (abs(tof32(1,i,iadc)).gt.10000.) tof32(1,i,iadc)= 10000.
         if (abs(tof32(2,i,iadc)).gt.10000.) tof32(2,i,iadc)= 10000.
      ENDDO

C----------------------------------------------------------------------

C------read  tracking routine
*     igoodevent = igoodevent+1
*     assigned input  parameters for track routine
*     1) Z-coordinates where the trajectory is evaluated
      do itof=1,NPTOF
         ZIN(itof) = ZTOF(itof)
      enddo
*     2) track status vector
C     COPY THE ALFA VECTOR FROM AL_PP TO AL_P FOR THE TRACK "T"
      do i=1,5
         AL_P(i) = al_pp(i)
      enddo

c      write(*,*) AL_P
c      write(*,*) 'Rigidity ',(1./AL_P(5))

      if (al_p(5).eq.0.) THEN
c         PRINT *,' TOF - WARNING F77: track with R = 0, discarded'
         GOTO 969
      ENDIF
*     -------- *** tracking routine *** --------
      IFAIL = 0
C     call DOTRACK(NPTOF,ZIN,XOUT,YOUT,AL_P,IFAIL)
      call DOTRACK2(NPTOF,ZIN,XOUT,YOUT,THXOUT,THYOUT,TLOUT,AL_P,IFAIL)

C     write(*,*) (TLOUT(i),i=1,6)

      if(IFAIL.ne.0)then
c         print *,' TOF - WARNING F77: tracking failed '
         goto 969
      endif
*     ------------------------------------------

 969  continue

C---  Fill xtr_tof  and ytr_tof: positions from tracker at ToF layers
      do j=1,6
      xtr_tof(j) = XOUT(j)
      ytr_tof(j) = YOUT(j)
c      write(*,*) XOUT(j),YOUT(j)
      enddo


C---  convert  angles to radian
      do j=1,6
      THXOUT(j) = 3.1415927*THXOUT(j)/180.
      THYOUT(j) = 3.1415927*THYOUT(j)/180.
      enddo

      do j=1,6
c      write (*,*) j,THXOUT(j),THYOUT(j)
      enddo


C----------------------------------------------------------------------
C------------------  set ADC & TDC flag = 0    ------------------------
C----------------------------------------------------------------------

      do j=1,8
      if (adc(ch11a(j),hb11a(j)).LT.4096)adcflagtof(ch11a(j),hb11a(j))=0
      if (adc(ch11b(j),hb11b(j)).LT.4096)adcflagtof(ch11b(j),hb11b(j))=0
      if (tdc(ch11a(j),hb11a(j)).LT.4096)tdcflagtof(ch11a(j),hb11a(j))=0
      if (tdc(ch11b(j),hb11b(j)).LT.4096)tdcflagtof(ch11b(j),hb11b(j))=0
      enddo
      do j=1,6
      if (adc(ch12a(j),hb12a(j)).LT.4096)adcflagtof(ch12a(j),hb12a(j))=0
      if (adc(ch12b(j),hb12b(j)).LT.4096)adcflagtof(ch12b(j),hb12b(j))=0
      if (tdc(ch12a(j),hb12a(j)).LT.4096)tdcflagtof(ch12a(j),hb12a(j))=0
      if (tdc(ch12b(j),hb12b(j)).LT.4096)tdcflagtof(ch12b(j),hb12b(j))=0
      enddo
      do j=1,2
      if (adc(ch21a(j),hb21a(j)).LT.4096)adcflagtof(ch21a(j),hb21a(j))=0
      if (adc(ch21b(j),hb21b(j)).LT.4096)adcflagtof(ch21b(j),hb21b(j))=0
      if (tdc(ch21a(j),hb21a(j)).LT.4096)tdcflagtof(ch21a(j),hb21a(j))=0
      if (tdc(ch21b(j),hb21b(j)).LT.4096)tdcflagtof(ch21b(j),hb21b(j))=0
      enddo
      do j=1,2
      if (adc(ch22a(j),hb22a(j)).LT.4096)adcflagtof(ch22a(j),hb22a(j))=0
      if (adc(ch22b(j),hb22b(j)).LT.4096)adcflagtof(ch22b(j),hb22b(j))=0
      if (tdc(ch22a(j),hb22a(j)).LT.4096)tdcflagtof(ch22a(j),hb22a(j))=0
      if (tdc(ch22b(j),hb22b(j)).LT.4096)tdcflagtof(ch22b(j),hb22b(j))=0
      enddo
      do j=1,3
      if (adc(ch31a(j),hb31a(j)).LT.4096)adcflagtof(ch31a(j),hb31a(j))=0
      if (adc(ch31b(j),hb31b(j)).LT.4096)adcflagtof(ch31b(j),hb31b(j))=0
      if (tdc(ch31a(j),hb31a(j)).LT.4096)tdcflagtof(ch31a(j),hb31a(j))=0
      if (tdc(ch31b(j),hb31b(j)).LT.4096)tdcflagtof(ch31b(j),hb31b(j))=0
      enddo
      do j=1,3
      if (adc(ch32a(j),hb32a(j)).LT.4096)adcflagtof(ch32a(j),hb32a(j))=0
      if (adc(ch32b(j),hb32b(j)).LT.4096)adcflagtof(ch32b(j),hb32b(j))=0
      if (tdc(ch32a(j),hb32a(j)).LT.4096)tdcflagtof(ch32a(j),hb32a(j))=0
      if (tdc(ch32b(j),hb32b(j)).LT.4096)tdcflagtof(ch32b(j),hb32b(j))=0
      enddo


C----------------------------------------------------------------
C---------- Check PMTs 10 and 35 for strange TDC values----------
C----------------------------------------------------------------

C---- S116A TDC=819
       if (tof11(1,6,1).EQ.819) then
             tof11(1,6,1) = 4095
             tdcflagtof(ch11a(6),hb11a(6))=2
       endif

C---- S222B TDC=819
       if (tof22(2,2,1).EQ.819) then
             tof22(2,2,1) = 4095
             tdcflagtof(ch22b(2),hb22b(2))=2
       endif

C-------------------------------------------------------------
C-------check which paddle penetrated the track  -----------
C-------------------------------------------------------------
c     middle y (or x) position of the upper and middle ToF-Paddle
c     DATA tof11_x/ -17.85,-12.75,-7.65,-2.55,2.55,7.65,12.75,17.85/
c     DATA tof12_y/ -13.75,-8.25,-2.75,2.75,8.25,13.75/
c     DATA tof21_y/  3.75,-3.75/  ! paddles in different order
c     DATA tof22_x/ -4.5,4.5/
c     DATA tof31_x/ -6.0,0.,6.0/
c     DATA tof32_y/ -5.0,0.0,5.0/
c
c     S11 8 paddles  33.0 x 5.1 cm
c     S12 6 paddles  40.8 x 5.5 cm
c     S21 2 paddles  18.0 x 7.5 cm
c     S22 2 paddles  15.0 x 9.0 cm
c     S31 3 paddles  15.0 x 6.0 cm
c     S32 3 paddles  18.0 x 5.0 cm


C--------------S11 --------------------------------------

      tof11_i = none_find

      yl = -33.0/2.
      yh =  33.0/2.

      if ((yout(1).gt.yl).and.(yout(1).lt.yh)) then
         do i=1,8
            xl = tof11_x(i) - 5.1/2.
            xh = tof11_x(i) + 5.1/2.
            if ((xout(1).gt.xl).and.(xout(1).le.xh)) then
               tof11_i=i
            endif
         enddo
      endif

C--------------S12 --------------------------------------

      tof12_i = none_find

      xl = -40.8/2.
      xh =  40.8/2.

      if ((xout(2).gt.xl).and.(xout(2).lt.xh)) then
         do i=1,6
            yl = tof12_y(i) - 5.5/2.
            yh = tof12_y(i) + 5.5/2.
            if ((yout(2).gt.yl).and.(yout(2).le.yh)) then
               tof12_i=i
            endif
         enddo
      endif

C--------------S21 --------------------------------------

      tof21_i = none_find

      xl = -18./2.
      xh =  18./2.

      if ((xout(3).gt.xl).and.(xout(3).lt.xh)) then
         do i=1,2
            yl = tof21_y(i) - 7.5/2.
            yh = tof21_y(i) + 7.5/2.
            if ((yout(3).gt.yl).and.(yout(3).le.yh)) then
               tof21_i=i
            endif
         enddo
      endif

C--------------S22 --------------------------------------

      tof22_i = none_find

      yl = -15./2.
      yh =  15./2.

      if ((yout(4).gt.yl).and.(yout(4).lt.yh)) then
         do i=1,2
            xl = tof22_x(i) - 9.0/2.
            xh = tof22_x(i) + 9.0/2.
            if ((xout(4).gt.xl).and.(xout(4).le.xh)) then
               tof22_i=i
            endif
         enddo
      endif

C--------------S31 --------------------------------------

      tof31_i = none_find

      yl = -15.0/2.
      yh =  15.0/2.

      if ((yout(5).gt.yl).and.(yout(5).lt.yh)) then
         do i=1,3
            xl = tof31_x(i) - 6.0/2.
            xh = tof31_x(i) + 6.0/2.
            if ((xout(5).gt.xl).and.(xout(5).le.xh)) then
               tof31_i=i
            endif
         enddo
      endif

C--------------S32 --------------------------------------

      tof32_i = none_find

      xl = -18.0/2.
      xh =  18.0/2.

      if ((xout(6).gt.xl).and.(xout(6).lt.xh)) then
         do i=1,3
            yl = tof32_y(i) - 5.0/2.
            yh = tof32_y(i) + 5.0/2.
            if ((yout(6).gt.yl).and.(yout(6).le.yh)) then
               tof32_i=i
            endif
         enddo
      endif


       hitvec(1)=tof11_i
       hitvec(2)=tof12_i
       hitvec(3)=tof21_i
       hitvec(4)=tof22_i
       hitvec(5)=tof31_i
       hitvec(6)=tof32_i


C----------------------------------------------------------------------
C--- check charge:
C--- if Z=2 we should use the attenuation curve for helium to
C--- fill the artificail ADC values and NOT divide by "hepratio"
C--- if Z>2 we should do a correction to
C--- the k1 constants in the beta calculation
C----------------------------------------------------------------------

       theta=0.   
       dist = ZTOF(1) - ZTOF(5)
       dl = 0.
       DO I=1,5
         dl = dl + TLOUT(i)
       ENDDO
       F = dl/dist
       theta = acos(1/F)

       iz = int(check_charge(theta,hitvec))
c       write(*,*) 'in toftrk',iz


C-------------------------------  new  ---------------------------
C--  calculate track position in paddle using timing difference
C--  this calculation is preliminary and uses some standard
C--  calibration values, but we need to find a rough position to
C--  be able to calculate artificial ADC values (needed for the
C--  timewalk...
C------------------------------------------------------------------

       do i=1,3
         xtofpre(i)=100.
         ytofpre(i)=100.
       enddo

C-----------------------------S1 --------------------------------

      IF (tof11_i.GT.none_find) THEN
      IF ((tof11(1,tof11_i,itdc).LT.2000).AND.
     +                             (tof11(2,tof11_i,itdc).LT.2000))
     +    ytofpre(1)  = ((tof11(1,tof11_i,itdc)-tof11(2,tof11_i,itdc))/2.
     +   -y_coor_lin11c(tof11_i,offset))/y_coor_lin11c(tof11_i,slope)
      endif

      IF (tof12_i.GT.none_find) THEN
      IF ((tof12(1,tof12_i,itdc).LT.2000).AND.
     +                             (tof12(2,tof12_i,itdc).LT.2000))
     +    xtofpre(1)  = ((tof12(1,tof12_i,itdc)-tof12(2,tof12_i,itdc))/2.
     +   -x_coor_lin12c(tof12_i,offset))/x_coor_lin12c(tof12_i,slope)
      endif


C-----------------------------S2 --------------------------------

      IF (tof21_i.GT.none_find) THEN
      IF ((tof21(1,tof21_i,itdc).LT.2000).AND.
     +                             (tof21(2,tof21_i,itdc).LT.2000))
     +    xtofpre(2) = ((tof21(1,tof21_i,itdc)-tof21(2,tof21_i,itdc))/2.
     +    -x_coor_lin21c(tof21_i,offset))/x_coor_lin21c(tof21_i,slope)
      endif

      IF (tof22_i.GT.none_find) THEN
      IF ((tof22(1,tof22_i,itdc).LT.2000).AND.
     +                             (tof22(2,tof22_i,itdc).LT.2000))
     +    ytofpre(2) = ((tof22(1,tof22_i,itdc)-tof22(2,tof22_i,itdc))/2.
     +    -y_coor_lin22c(tof22_i,offset))/y_coor_lin22c(tof22_i,slope)
      endif


C-----------------------------S3 --------------------------------

      IF (tof31_i.GT.none_find) THEN
      IF ((tof31(1,tof31_i,itdc).LT.2000).AND.
     +                             (tof31(2,tof31_i,itdc).LT.2000))
     +    ytofpre(3)  = ((tof31(1,tof31_i,itdc)-tof31(2,tof31_i,itdc))/2.
     +    -y_coor_lin31c(tof31_i,offset))/y_coor_lin31c(tof31_i,slope)
      endif

      IF (tof32_i.GT.none_find) THEN
      IF ((tof32(1,tof32_i,itdc).LT.2000).AND.
     +                             (tof32(2,tof32_i,itdc).LT.2000))
     +    xtofpre(3)  = ((tof32(1,tof32_i,itdc)-tof32(2,tof32_i,itdc))/2.
     +    -x_coor_lin32c(tof32_i,offset))/x_coor_lin32c(tof32_i,slope)
      endif


C--  restrict TDC measurements to physical paddle dimensions +/- 10 cm

        if (abs(xtofpre(1)).gt.31.)  xtofpre(1)=100.
        if (abs(xtofpre(2)).gt.19.)  xtofpre(2)=100.
        if (abs(xtofpre(3)).gt.19.)  xtofpre(3)=100.

        if (abs(ytofpre(1)).gt.26.)  ytofpre(1)=100.
        if (abs(ytofpre(2)).gt.18.)  ytofpre(2)=100.
        if (abs(ytofpre(3)).gt.18.)  ytofpre(3)=100.

C--------------------------------------------------------------------
C---- if paddle hit: if we have TDC value but no ADC, create ADC value
C---- use the "pre" position if possible, since this gives better time
C---- resolution ... october 2008
C--------------------------------------------------------------------
c     middle y (or x) position of the upper and middle ToF-Paddle
c     DATA tof11_x/ -17.85,-12.75,-7.65,-2.55,2.55,7.65,12.75,17.85/
c     DATA tof12_y/ -13.75,-8.25,-2.75,2.75,8.25,13.75/
c     DATA tof21_y/  3.75,-3.75/  ! paddles in different order
c     DATA tof22_x/ -4.5,4.5/
c     DATA tof31_x/ -6.0,0.,6.0/
c     DATA tof32_y/ -5.0,0.0,5.0/

C----------------------------S1 -------------------------------------

c     yhelp=yout(1)
      yhelp = ytofpre(1)
      if (yhelp.eq.100) yhelp=yout(1)

      IF (tof11_i.GT.none_find.AND.abs(yhelp).lt.100) THEN
         i = tof11_i
         if ((tdc(ch11a(i),hb11a(i)).lt.4095).AND.
     &       (adc(ch11a(i),hb11a(i)).eq.4095)) then
            phi   = atan(tan(THYOUT(1))/tan(THXOUT(1)))
            theta = atan(tan(THXOUT(1))/cos(phi))
            xkorr = atten(left,11,i,yhelp)
            if (iz.le.1) xkorr=xkorr/hepratio
            tof11(left,i,iadc)=xkorr/cos(theta)
            adcflag(ch11a(i),hb11a(i)) = 1
         endif
         if ((tdc(ch11b(i),hb11b(i)).lt.4095).AND.
     &       (adc(ch11b(i),hb11b(i)).eq.4095)) then
            phi   = atan(tan(THYOUT(1))/tan(THXOUT(1)))
            theta = atan(tan(THXOUT(1))/cos(phi))
            xkorr = atten(right,11,i,yhelp)
            if (iz.le.1) xkorr=xkorr/hepratio
            tof11(right,i,iadc)=xkorr/cos(theta)
            adcflag(ch11b(i),hb11b(i)) = 1
         endif
      ENDIF

c      xhelp=xout(2)
      xhelp = xtofpre(1)
      if (xhelp.eq.100) xhelp=xout(2)

      IF (tof12_i.GT.none_find.AND.abs(xhelp).lt.100) THEN
         i = tof12_i
         if ((tdc(ch12a(i),hb12a(i)).lt.4095).AND.
     &       (adc(ch12a(i),hb12a(i)).eq.4095)) then
            phi   = atan(tan(THYOUT(2))/tan(THXOUT(2)))
            theta = atan(tan(THXOUT(2))/cos(phi))
c            xkorr=adcx12(left,i,1)*exp(-xhelp/adcx12(left,i,2))
            xkorr = atten(left,12,i,xhelp)
            if (iz.le.1) xkorr=xkorr/hepratio
            tof12(left,i,iadc) = xkorr/cos(theta)
            adcflag(ch12a(i),hb12a(i)) = 1
         endif
         if ((tdc(ch12b(i),hb12b(i)).lt.4095).AND.
     &       (adc(ch12b(i),hb12b(i)).eq.4095)) then
            phi   = atan(tan(THYOUT(2))/tan(THXOUT(2)))
            theta = atan(tan(THXOUT(2))/cos(phi))
c            xkorr=adcx12(right,i,1)*exp(xhelp/adcx12(right,i,2))
            xkorr = atten(right,12,i,xhelp)
            if (iz.le.1) xkorr=xkorr/hepratio
            tof12(right,i,iadc) = xkorr/cos(theta)
            adcflag(ch12b(i),hb12b(i)) = 1
         endif
      ENDIF

C-----------------------------S2 --------------------------------

c      xhelp=xout(3)
      xhelp = xtofpre(2)
      if (xhelp.eq.100) xhelp=xout(3)

      IF (tof21_i.GT.none_find.AND.abs(xhelp).lt.100) THEN
         i = tof21_i
         if ((tdc(ch21a(i),hb21a(i)).lt.4095).AND.
     &       (adc(ch21a(i),hb21a(i)).eq.4095)) then
            phi   = atan(tan(THYOUT(3))/tan(THXOUT(3)))
            theta = atan(tan(THXOUT(3))/cos(phi))
c            xkorr=adcx21(left,i,1)*exp(-xhelp/adcx21(left,i,2))
            xkorr = atten(left,21,i,xhelp)
            if (iz.le.1) xkorr=xkorr/hepratio
            tof21(left,i,iadc) = xkorr/cos(theta)
            adcflag(ch21a(i),hb21a(i)) = 1
         endif
         if ((tdc(ch21b(i),hb21b(i)).lt.4095).AND.
     &       (adc(ch21b(i),hb21b(i)).eq.4095)) then
            phi   = atan(tan(THYOUT(3))/tan(THXOUT(3)))
            theta = atan(tan(THXOUT(3))/cos(phi))
c            xkorr=adcx21(right,i,1)*exp(xhelp/adcx21(right,i,2))
            xkorr = atten(right,21,i,xhelp)
            if (iz.le.1) xkorr=xkorr/hepratio
            tof21(right,i,iadc) = xkorr/cos(theta)
            adcflag(ch21b(i),hb21b(i)) = 1
         endif
      ENDIF


c      yhelp=yout(4)
      yhelp = ytofpre(2)
      if (yhelp.eq.100) yhelp=yout(4)

      IF (tof22_i.GT.none_find.AND.abs(yhelp).lt.100) THEN
         i = tof22_i
         if ((tdc(ch22a(i),hb22a(i)).lt.4095).AND.
     &       (adc(ch22a(i),hb22a(i)).eq.4095)) then
            phi   = atan(tan(THYOUT(4))/tan(THXOUT(4)))
            theta = atan(tan(THXOUT(4))/cos(phi))
c            xkorr=adcx22(left,i,1)*exp(-yhelp/adcx22(left,i,2))
            xkorr = atten(left,22,i,yhelp)
            if (iz.le.1) xkorr=xkorr/hepratio
            tof22(left,i,iadc) = xkorr/cos(theta)
            adcflag(ch22a(i),hb22a(i)) = 1
         endif
         if ((tdc(ch22b(i),hb22b(i)).lt.4095).AND.
     &       (adc(ch22b(i),hb22b(i)).eq.4095)) then
            phi   = atan(tan(THYOUT(4))/tan(THXOUT(4)))
            theta = atan(tan(THXOUT(4))/cos(phi))
c            xkorr=adcx22(right,i,1)*exp(yhelp/adcx22(right,i,2))
            xkorr = atten(right,22,i,yhelp)
            if (iz.le.1) xkorr=xkorr/hepratio
            tof22(right,i,iadc) = xkorr/cos(theta)
            adcflag(ch22b(i),hb22b(i)) = 1
         endif
      ENDIF

C-----------------------------S3 --------------------------------

c      yhelp=yout(5)
      yhelp = ytofpre(3)
      if (yhelp.eq.100) yhelp=yout(5)

      IF (tof31_i.GT.none_find.AND.abs(yhelp).lt.100) THEN
         i = tof31_i
         if ((tdc(ch31a(i),hb31a(i)).lt.4095).AND.
     &       (adc(ch31a(i),hb31a(i)).eq.4095)) then
            phi   = atan(tan(THYOUT(5))/tan(THXOUT(5)))
            theta = atan(tan(THXOUT(5))/cos(phi))
c            xkorr=adcx31(left,i,1)*exp(-yhelp/adcx31(left,i,2))
            xkorr = atten(left,31,i,yhelp)
            if (iz.le.1) xkorr=xkorr/hepratio
            tof31(left,i,iadc) = xkorr/cos(theta)
            adcflag(ch31a(i),hb31a(i)) = 1
         endif
         if ((tdc(ch31b(i),hb31b(i)).lt.4095).AND.
     &       (adc(ch31b(i),hb31b(i)).eq.4095)) then
            phi   = atan(tan(THYOUT(5))/tan(THXOUT(5)))
            theta = atan(tan(THXOUT(5))/cos(phi))
c            xkorr=adcx31(right,i,1)*exp(yhelp/adcx31(right,i,2))
            xkorr = atten(right,31,i,yhelp)
            if (iz.le.1) xkorr=xkorr/hepratio
            tof31(right,i,iadc) = xkorr/cos(theta)
            adcflag(ch31b(i),hb31b(i)) = 1
         endif
      ENDIF


c      xhelp=xout(6)
      xhelp = xtofpre(3)
      if (xhelp.eq.100) xhelp=xout(6)

      IF (tof32_i.GT.none_find.AND.abs(xout(6)).lt.100) THEN
         i = tof32_i
         if ((tdc(ch32a(i),hb32a(i)).lt.4095).AND.
     &       (adc(ch32a(i),hb32a(i)).eq.4095)) then
            phi   = atan(tan(THYOUT(6))/tan(THXOUT(6)))
            theta = atan(tan(THXOUT(6))/cos(phi))
c            xkorr=adcx32(left,i,1)*exp(-xhelp/adcx32(left,i,2))
            xkorr = atten(left,32,i,xhelp)
            if (iz.le.1) xkorr=xkorr/hepratio
            tof32(left,i,iadc) = xkorr/cos(theta)
            adcflag(ch32a(i),hb32a(i)) = 1
         endif
         if ((tdc(ch32b(i),hb32b(i)).lt.4095).AND.
     &       (adc(ch32b(i),hb32b(i)).eq.4095)) then
            phi   = atan(tan(THYOUT(6))/tan(THXOUT(6)))
            theta = atan(tan(THXOUT(6))/cos(phi))
c            xkorr=adcx32(right,i,1)*exp(xhelp/adcx32(right,i,2))
            xkorr = atten(right,32,i,xhelp)
            if (iz.le.1) xkorr=xkorr/hepratio
            tof32(right,i,iadc) = xkorr/cos(theta)
            adcflag(ch32b(i),hb32b(i)) = 1
         endif
      ENDIF

C-------------------------------------------------------------------
C Now there is for each hitted paddle a TDC and ADC value, if the 
C TDC was < 4095. 
C There might be also TDC-ADC pairs in paddles not hitted
C Let's correct the raw TDC value with the time walk
C------------------------------------------------------------------- 
C--------------------Time walk correction  -------------------------
C-------------------------------------------------------------------

      DO i=1,8
         if ((tdc(ch11a(i),hb11a(i)).lt.4095).and.
     &             (tof11(left,i,iadc).lt.3786)) THEN
         xhelp = tw11(left,i)/(tof11(left,i,iadc)**0.5)
         tof11(left,i,itdc) = tof11(left,i,itdc) + xhelp
         tdc_c(ch11a(i),hb11a(i))=tof11(left,i,itdc)
                                              ENDIF 

         if ((tdc(ch11b(i),hb11b(i)).lt.4095).and.
     &             (tof11(right,i,iadc).lt.3786)) THEN
         xhelp = tw11(right,i)/(tof11(right,i,iadc)**0.5)
         tof11(right,i,itdc) = tof11(right,i,itdc) + xhelp
         tdc_c(ch11b(i),hb11b(i))=tof11(right,i,itdc)
                                             ENDIF
      ENDDO


      DO i=1,6
         if ((tdc(ch12a(i),hb12a(i)).lt.4095).and.
     &             (tof12(left,i,iadc).lt.3786)) THEN
         xhelp = tw12(left,i)/(tof12(left,i,iadc)**0.5)
         tof12(left,i,itdc) = tof12(left,i,itdc) + xhelp
         tdc_c(ch12a(i),hb12a(i))=tof12(left,i,itdc)
                                              ENDIF 

         if ((tdc(ch12b(i),hb12b(i)).lt.4095).and.
     &             (tof12(right,i,iadc).lt.3786)) THEN
         xhelp = tw12(right,i)/(tof12(right,i,iadc)**0.5)
         tof12(right,i,itdc) = tof12(right,i,itdc) + xhelp
         tdc_c(ch12b(i),hb12b(i))=tof12(right,i,itdc)
                                             ENDIF
      ENDDO

C----
      DO I=1,2 
         if ((tdc(ch21a(i),hb21a(i)).lt.4095).and.
     &             (tof21(left,i,iadc).lt.3786)) THEN
         xhelp = tw21(left,i)/(tof21(left,i,iadc)**0.5)
         tof21(left,i,itdc) = tof21(left,i,itdc) + xhelp
         tdc_c(ch21a(i),hb21a(i))=tof21(left,i,itdc)
                                              ENDIF 

         if ((tdc(ch21b(i),hb21b(i)).lt.4095).and.
     &             (tof21(right,i,iadc).lt.3786)) THEN
         xhelp = tw21(right,i)/(tof21(right,i,iadc)**0.5)
         tof21(right,i,itdc) = tof21(right,i,itdc) + xhelp
         tdc_c(ch21b(i),hb21b(i))=tof21(right,i,itdc)
                                             ENDIF
      ENDDO

      DO I=1,2 
         if ((tdc(ch22a(i),hb22a(i)).lt.4095).and.
     &             (tof22(left,i,iadc).lt.3786)) THEN
         xhelp = tw22(left,i)/(tof22(left,i,iadc)**0.5)
         tof22(left,i,itdc) = tof22(left,i,itdc) + xhelp
         tdc_c(ch22a(i),hb22a(i))=tof22(left,i,itdc)
                                              ENDIF 

         if ((tdc(ch22b(i),hb22b(i)).lt.4095).and.
     &             (tof22(right,i,iadc).lt.3786)) THEN
         xhelp = tw22(right,i)/(tof22(right,i,iadc)**0.5)
         tof22(right,i,itdc) = tof22(right,i,itdc) + xhelp
         tdc_c(ch22b(i),hb22b(i))=tof22(right,i,itdc)
                                             ENDIF
      ENDDO

C----
      DO I=1,3 
         if ((tdc(ch31a(i),hb31a(i)).lt.4095).and.
     &             (tof31(left,i,iadc).lt.3786)) THEN
         xhelp = tw31(left,i)/(tof31(left,i,iadc)**0.5)
         tof31(left,i,itdc) = tof31(left,i,itdc) + xhelp
         tdc_c(ch31a(i),hb31a(i))=tof31(left,i,itdc)
                                              ENDIF 

         if ((tdc(ch31b(i),hb31b(i)).lt.4095).and.
     &             (tof31(right,i,iadc).lt.3786)) THEN
         xhelp = tw31(right,i)/(tof31(right,i,iadc)**0.5)
         tof31(right,i,itdc) = tof31(right,i,itdc) + xhelp
         tdc_c(ch31b(i),hb31b(i))=tof31(right,i,itdc)
                                             ENDIF
      ENDDO

      DO I=1,3
         if ((tdc(ch32a(i),hb32a(i)).lt.4095).and.
     &             (tof32(left,i,iadc).lt.3786)) THEN
         xhelp = tw32(left,i)/(tof32(left,i,iadc)**0.5)
         tof32(left,i,itdc) = tof32(left,i,itdc) + xhelp
         tdc_c(ch32a(i),hb32a(i))=tof32(left,i,itdc)
                                              ENDIF 

         if ((tdc(ch32b(i),hb32b(i)).lt.4095).and.
     &             (tof32(right,i,iadc).lt.3786)) THEN
         xhelp = tw32(right,i)/(tof32(right,i,iadc)**0.5)
         tof32(right,i,itdc) = tof32(right,i,itdc) + xhelp
         tdc_c(ch32b(i),hb32b(i))=tof32(right,i,itdc)
                                             ENDIF
      ENDDO


C-----------------------------------------------------------------------
C--------------------Insert Artifical TDC Value  ---------------------
C     For each Paddle perform check:
C     if left paddle=4095  and right paddle OK => create TDC value left
C     if right paddle=4095  and left paddle OK => create TDC value right
C-----------------------------------------------------------------------

C-----------------------S11 -----------------------------------------

      IF (tof11_i.GT.none_find) THEN
         xpos = yout(1)
         i = tof11_i
        if ((tdc(ch11a(i),hb11a(i)).EQ.4095).AND.
     &          (tdc(ch11b(i),hb11b(i)).LT.4095)) THEN
          tof11(1,tof11_i,itdc) = tof11(2,tof11_i,itdc)
     &           + 2*(y_coor_lin11(tof11_i,offset)
     &           + xpos*y_coor_lin11(tof11_i,slope))
            tdcflag(ch11a(i),hb11a(i)) = 1
        ENDIF

        if ((tdc(ch11b(i),hb11b(i)).EQ.4095).AND.
     &          (tdc(ch11a(i),hb11a(i)).LT.4095)) THEN
            tof11(2,tof11_i,itdc) = tof11(1,tof11_i,itdc)
     &           - 2*(y_coor_lin11(tof11_i,offset)
     &           + xpos*y_coor_lin11(tof11_i,slope))
            tdcflag(ch11b(i),hb11b(i)) = 1
         ENDIF

      ENDIF

C-----------------------S12 -----------------------------------------

      IF (tof12_i.GT.none_find) THEN
         xpos = xout(2)
         i = tof12_i
        if ((tdc(ch12a(i),hb12a(i)).EQ.4095).AND.
     &          (tdc(ch12b(i),hb12b(i)).LT.4095)) THEN
            tof12(1,tof12_i,itdc) = tof12(2,tof12_i,itdc)
     &           + 2*(x_coor_lin12(tof12_i,offset)
     &           + xpos*x_coor_lin12(tof12_i,slope))
            tdcflag(ch12a(i),hb12a(i)) = 1
         ENDIF

        if ((tdc(ch12b(i),hb12b(i)).EQ.4095).AND.
     &          (tdc(ch12a(i),hb12a(i)).LT.4095)) THEN
            tof12(2,tof12_i,itdc) = tof12(1,tof12_i,itdc)
     &           - 2*(x_coor_lin12(tof12_i,offset)
     &           + xpos*x_coor_lin12(tof12_i,slope))
            tdcflag(ch12b(i),hb12b(i)) = 1
         ENDIF
      ENDIF

C-----------------------S21 -----------------------------------------

      IF (tof21_i.GT.none_find) THEN
         xpos = xout(3)
         i = tof21_i
        if ((tdc(ch21a(i),hb21a(i)).EQ.4095).AND.
     &          (tdc(ch21b(i),hb21b(i)).LT.4095)) THEN
            tof21(1,tof21_i,itdc) = tof21(2,tof21_i,itdc)
     &           + 2*(x_coor_lin21(tof21_i,offset)
     &           + xpos*x_coor_lin21(tof21_i,slope))
            tdcflag(ch21a(i),hb21a(i)) = 1
         ENDIF

        if ((tdc(ch21b(i),hb21b(i)).EQ.4095).AND.
     &          (tdc(ch21a(i),hb21a(i)).LT.4095)) THEN
            tof21(2,tof21_i,itdc) = tof21(1,tof21_i,itdc)
     &           - 2*(x_coor_lin21(tof21_i,offset)
     &           + xpos*x_coor_lin21(tof21_i,slope))
            tdcflag(ch21b(i),hb21b(i)) = 1
         ENDIF
      ENDIF

C-----------------------S22 -----------------------------------------

      IF (tof22_i.GT.none_find) THEN
         xpos = yout(4)
         i = tof22_i
        if ((tdc(ch22a(i),hb22a(i)).EQ.4095).AND.
     &          (tdc(ch22b(i),hb22b(i)).LT.4095)) THEN
            tof22(1,tof22_i,itdc) = tof22(2,tof22_i,itdc)
     &           + 2*(y_coor_lin22(tof22_i,offset)
     &           + xpos*y_coor_lin22(tof22_i,slope))
            tdcflag(ch22a(i),hb22a(i)) = 1
         ENDIF

        if ((tdc(ch22b(i),hb22b(i)).EQ.4095).AND.
     &          (tdc(ch22a(i),hb22a(i)).LT.4095)) THEN
            tof22(2,tof22_i,itdc) = tof22(1,tof22_i,itdc)
     &           - 2*(y_coor_lin22(tof22_i,offset)
     &           + xpos*y_coor_lin22(tof22_i,slope))
            tdcflag(ch22b(i),hb22b(i)) = 1
         ENDIF
      ENDIF

C-----------------------S31 -----------------------------------------

      IF (tof31_i.GT.none_find) THEN
         xpos = yout(5)
         i = tof31_i
        if ((tdc(ch31a(i),hb31a(i)).EQ.4095).AND.
     &          (tdc(ch31b(i),hb31b(i)).LT.4095)) THEN
            tof31(1,tof31_i,itdc) = tof31(2,tof31_i,itdc)
     &           + 2*(y_coor_lin31(tof31_i,offset)
     &           + xpos*y_coor_lin31(tof31_i,slope))
            tdcflag(ch31a(i),hb31a(i)) = 1
         ENDIF

        if ((tdc(ch31b(i),hb31b(i)).EQ.4095).AND.
     &          (tdc(ch31a(i),hb31a(i)).LT.4095)) THEN
            tof31(2,tof31_i,itdc) = tof31(1,tof31_i,itdc)
     &           - 2*(y_coor_lin31(tof31_i,offset)
     &           + xpos*y_coor_lin31(tof31_i,slope))
            tdcflag(ch31b(i),hb31b(i)) = 1
         ENDIF
      ENDIF

C-----------------------S32 -----------------------------------------

      IF (tof32_i.GT.none_find) THEN
         xpos = xout(6)
         i = tof32_i
        if ((tdc(ch32a(i),hb32a(i)).EQ.4095).AND.
     &          (tdc(ch32b(i),hb32b(i)).LT.4095)) THEN
            tof32(1,tof32_i,itdc) = tof32(2,tof32_i,itdc)
     &           + 2*(x_coor_lin32(tof32_i,offset)
     &           + xpos*x_coor_lin32(tof32_i,slope))
            tdcflag(ch32a(i),hb32a(i)) = 1
         ENDIF

        if ((tdc(ch32b(i),hb32b(i)).EQ.4095).AND.
     &          (tdc(ch32a(i),hb32a(i)).LT.4095)) THEN
            tof32(2,tof32_i,itdc) = tof32(1,tof32_i,itdc)
     &           - 2*(x_coor_lin32(tof32_i,offset)
     &           + xpos*x_coor_lin32(tof32_i,slope))
            tdcflag(ch32b(i),hb32b(i)) = 1
         ENDIF
      ENDIF


C------------------------------------------------------------------
C---  calculate track position in paddle using timing difference
C------------------------------------------------------------------

      do i=1,3
         xtofpos(i)=100.
         ytofpos(i)=100.
      enddo

C-----------------------------S1 --------------------------------

      IF (tof11_i.GT.none_find) THEN
         IF ((tof11(1,tof11_i,itdc).NE.4095).AND.
     &       (tof11(2,tof11_i,itdc).NE.4095)) THEN
         ytofpos(1)  = ((tof11(1,tof11_i,itdc)-tof11(2,tof11_i,itdc))/2.
     +        -y_coor_lin11(tof11_i,offset))/y_coor_lin11(tof11_i,slope)
        if (abs(ytofpos(1)).gt.26.)  ytofpos(1)=101.
      i=tof11_i
      endif
      endif

      IF (tof12_i.GT.none_find) THEN
         IF ((tof12(1,tof12_i,itdc).NE.4095).AND.
     &       (tof12(2,tof12_i,itdc).NE.4095)) THEN
         xtofpos(1)  = ((tof12(1,tof12_i,itdc)-tof12(2,tof12_i,itdc))/2.
     +        -x_coor_lin12(tof12_i,offset))/x_coor_lin12(tof12_i,slope)
        if (abs(xtofpos(1)).gt.31.)  xtofpos(1)=101.
      i=tof12_i
      endif
      endif

C-----------------------------S2 --------------------------------

      IF (tof21_i.GT.none_find) THEN
         IF ((tof21(1,tof21_i,itdc).NE.4095).AND.
     &       (tof21(2,tof21_i,itdc).NE.4095)) THEN
         xtofpos(2) = ((tof21(1,tof21_i,itdc)-tof21(2,tof21_i,itdc))/2.
     +        -x_coor_lin21(tof21_i,offset))/x_coor_lin21(tof21_i,slope)
        if (abs(xtofpos(2)).gt.19.)  xtofpos(2)=101.
      i=tof21_i
      endif
      endif

      IF (tof22_i.GT.none_find) THEN
         IF ((tof22(1,tof22_i,itdc).NE.4095).AND.
     &       (tof22(2,tof22_i,itdc).NE.4095)) THEN
         ytofpos(2) = ((tof22(1,tof22_i,itdc)-tof22(2,tof22_i,itdc))/2.
     +        -y_coor_lin22(tof22_i,offset))/y_coor_lin22(tof22_i,slope)
        if (abs(ytofpos(2)).gt.18.)  ytofpos(2)=101.
      i=tof22_i
      endif
      endif

C-----------------------------S3 --------------------------------

      IF (tof31_i.GT.none_find) THEN
         IF ((tof31(1,tof31_i,itdc).NE.4095).AND.
     &       (tof31(2,tof31_i,itdc).NE.4095)) THEN
         ytofpos(3)  = ((tof31(1,tof31_i,itdc)-tof31(2,tof31_i,itdc))/2.
     +        -y_coor_lin31(tof31_i,offset))/y_coor_lin31(tof31_i,slope)
        if (abs(ytofpos(3)).gt.18.)  ytofpos(3)=101.
      i=tof31_i
      endif
      endif

      IF (tof32_i.GT.none_find) THEN
         IF ((tof32(1,tof32_i,itdc).NE.4095).AND.
     &       (tof32(2,tof32_i,itdc).NE.4095)) THEN
         xtofpos(3)  = ((tof32(1,tof32_i,itdc)-tof32(2,tof32_i,itdc))/2.
     +        -x_coor_lin32(tof32_i,offset))/x_coor_lin32(tof32_i,slope)
        if (abs(xtofpos(3)).gt.19.)  xtofpos(3)=101.
      i=tof32_i
      endif
      endif

c      do i=1,3
c         if (abs(xtofpos(i)).gt.100.) then
c            xtofpos(i)=101.
c         endif
c         if (abs(ytofpos(i)).gt.100.) then
c            ytofpos(i)=101.
c         endif
c      enddo


C--------------------------------------------------------------------
C-------------------Corrections on ADC-data -------------------------
C-----------------angle and ADC(x) correction -----------------------
C----------------   moved to the new dEdx routine -------------------

C--------------------------------------------------------------------
C----------------------calculate Beta  ------------------------------
C--------------------------------------------------------------------
C---------------------difference of sums  ---------------------------
C
C     DS = (t1+t2) - t3+t4)
C     DS = c1 + c2/beta*cos(theta)
C     c2 = 2d/c   gives c2 = 2d/(c*TDCresolution)  TDC=50ps/channel
C     =>  c2 =  ca.60 for  0.45 m    c2 = ca.109  for 0.81 m
C     since TDC resolution varies slightly c2 has to be calibrated
C     instead of cos(theta) use factor F:
C     F =  pathlength/d
C     => beta = c2*F/(DS-c1))

      dist = ZTOF(1) - ZTOF(5)
      dl = 0.
      DO I=1,5
         dl = dl + TLOUT(i)
      ENDDO
      F = dl/dist

C     S11 - S31
C      IF (tof11_i.GT.none_find.AND.tof31_i.GT.none_find) THEN
      IF ((tof11_i.GT.none_find).AND.(tof31_i.GT.none_find).AND.
     &    (ytofpos(1).NE.101.).AND.(ytofpos(3).NE.101.)) THEN
         t1 = tof11(1,tof11_i,itdc)
         t2 = tof11(2,tof11_i,itdc)
         t3 = tof31(1,tof31_i,itdc)
         t4 = tof31(2,tof31_i,itdc)
         IF ((t1.lt.4095).and.(t2.lt.4095).and.
     &        (t3.lt.4095).and.(t4.lt.4095)) THEN
            xhelp1 = tof11(1,tof11_i,itdc)+tof11(2,tof11_i,itdc)
            xhelp2 = tof31(1,tof31_i,itdc)+tof31(2,tof31_i,itdc)
            ds = xhelp1-xhelp2
            ihelp=(tof11_i-1)*3+tof31_i
            c1 = k_S11S31(1,ihelp)
            if (iz.gt.2) c1 = c1 + k1corrA1
            c2 = k_S11S31(2,ihelp)
            beta_a(1) = c2*F/(ds-c1)
c       write(*,*) 'S11-S31 ',c1,c2,F
c       write(*,*) 'S11-S31 ',xhelp1,xhelp2, beta_a(1)
C-------ToF Mask - S11 - S31

            tofmask(ch11a(tof11_i),hb11a(tof11_i)) =
     $           tofmask(ch11a(tof11_i),hb11a(tof11_i)) + 1
            tofmask(ch11b(tof11_i),hb11b(tof11_i)) =
     $           tofmask(ch11b(tof11_i),hb11b(tof11_i)) + 1

            tofmask(ch31a(tof31_i),hb31a(tof31_i)) =
     $           tofmask(ch31a(tof31_i),hb31a(tof31_i)) + 1
            tofmask(ch31b(tof31_i),hb31b(tof31_i)) =
     $           tofmask(ch31b(tof31_i),hb31b(tof31_i)) + 1

         ENDIF
      ENDIF

      dist = ZTOF(1) - ZTOF(6)
      dl = 0.
      DO I=1,6
         dl = dl + TLOUT(i)
      ENDDO
      F = dl/dist

C     S11 - S32
C      IF (tof11_i.GT.none_find.AND.tof32_i.GT.none_find) THEN
       IF ((tof11_i.GT.none_find).AND.(tof32_i.GT.none_find).AND.
     &    (ytofpos(1).NE.101.).AND.(xtofpos(3).NE.101.)) THEN
         t1 = tof11(1,tof11_i,itdc)
         t2 = tof11(2,tof11_i,itdc)
         t3 = tof32(1,tof32_i,itdc)
         t4 = tof32(2,tof32_i,itdc)
         IF ((t1.lt.4095).and.(t2.lt.4095).and.
     &        (t3.lt.4095).and.(t4.lt.4095)) THEN
            xhelp1 = tof11(1,tof11_i,itdc)+tof11(2,tof11_i,itdc)
            xhelp2 = tof32(1,tof32_i,itdc)+tof32(2,tof32_i,itdc)
            ds = xhelp1-xhelp2
            ihelp=(tof11_i-1)*3+tof32_i
            c1 = k_S11S32(1,ihelp)
            if (iz.gt.2) c1 = c1 + k1corrA1
            c2 = k_S11S32(2,ihelp)
            beta_a(2) = c2*F/(ds-c1)
C     write(*,*) 'S11-S32 ',xhelp1,xhelp2, beta_a(2)

C-------ToF Mask - S11 - S32

            tofmask(ch11a(tof11_i),hb11a(tof11_i)) =
     $           tofmask(ch11a(tof11_i),hb11a(tof11_i)) + 1
            tofmask(ch11b(tof11_i),hb11b(tof11_i)) =
     $           tofmask(ch11b(tof11_i),hb11b(tof11_i)) + 1

            tofmask(ch32a(tof32_i),hb32a(tof32_i)) =
     $           tofmask(ch32a(tof32_i),hb32a(tof32_i)) + 1
            tofmask(ch32b(tof32_i),hb32b(tof32_i)) =
     $           tofmask(ch32b(tof32_i),hb32b(tof32_i)) + 1

C-------

         ENDIF
      ENDIF

C     S12 - S31
      dist = ZTOF(2) - ZTOF(5)
      dl = 0.
      DO I=2,5
         dl = dl + TLOUT(i)
      ENDDO
      F = dl/dist

C      IF (tof12_i.GT.none_find.AND.tof31_i.GT.none_find) THEN
       IF ((tof12_i.GT.none_find).AND.(tof31_i.GT.none_find).AND.
     &    (xtofpos(1).NE.101.).AND.(ytofpos(3).NE.101.)) THEN
         t1 = tof12(1,tof12_i,itdc)
         t2 = tof12(2,tof12_i,itdc)
         t3 = tof31(1,tof31_i,itdc)
         t4 = tof31(2,tof31_i,itdc)
         IF ((t1.lt.4095).and.(t2.lt.4095).and.
     &        (t3.lt.4095).and.(t4.lt.4095)) THEN
            xhelp1 = tof12(1,tof12_i,itdc)+tof12(2,tof12_i,itdc)
            xhelp2 = tof31(1,tof31_i,itdc)+tof31(2,tof31_i,itdc)
            ds = xhelp1-xhelp2
            ihelp=(tof12_i-1)*3+tof31_i
            c1 = k_S12S31(1,ihelp)
            if (iz.gt.2) c1 = c1 + k1corrA1
            c2 = k_S12S31(2,ihelp)
            beta_a(3) = c2*F/(ds-c1)
C     write(*,*) 'S12-S31 ',xhelp1,xhelp2, beta_a(3)

C-------ToF Mask - S12 - S31

            tofmask(ch12a(tof12_i),hb12a(tof12_i)) =
     $           tofmask(ch12a(tof12_i),hb12a(tof12_i)) + 1
            tofmask(ch12b(tof12_i),hb12b(tof12_i)) =
     $           tofmask(ch12b(tof12_i),hb12b(tof12_i)) + 1

            tofmask(ch31a(tof31_i),hb31a(tof31_i)) =
     $           tofmask(ch31a(tof31_i),hb31a(tof31_i)) + 1
            tofmask(ch31b(tof31_i),hb31b(tof31_i)) =
     $           tofmask(ch31b(tof31_i),hb31b(tof31_i)) + 1

C-------

         ENDIF
      ENDIF

C     S12 - S32

      dist = ZTOF(2) - ZTOF(6)
      dl = 0.
      DO I=2,6
         dl = dl + TLOUT(i)
      ENDDO
      F = dl/dist

C      IF (tof12_i.GT.none_find.AND.tof32_i.GT.none_find) THEN
       IF ((tof12_i.GT.none_find).AND.(tof32_i.GT.none_find).AND.
     &    (xtofpos(1).NE.101.).AND.(xtofpos(3).NE.101.)) THEN
         t1 = tof12(1,tof12_i,itdc)
         t2 = tof12(2,tof12_i,itdc)
         t3 = tof32(1,tof32_i,itdc)
         t4 = tof32(2,tof32_i,itdc)
         IF ((t1.lt.4095).and.(t2.lt.4095).and.
     &        (t3.lt.4095).and.(t4.lt.4095)) THEN
            xhelp1 = tof12(1,tof12_i,itdc)+tof12(2,tof12_i,itdc)
            xhelp2 = tof32(1,tof32_i,itdc)+tof32(2,tof32_i,itdc)
            ds = xhelp1-xhelp2
            ihelp=(tof12_i-1)*3+tof32_i
            c1 = k_S12S32(1,ihelp)
            if (iz.gt.2) c1 = c1 + k1corrA1
            c2 = k_S12S32(2,ihelp)
            beta_a(4) = c2*F/(ds-c1)
C     write(*,*) 'S12-S32 ',xhelp1,xhelp2, beta_a(4)

C-------ToF Mask - S12 - S32

            tofmask(ch12a(tof12_i),hb12a(tof12_i)) =
     $           tofmask(ch12a(tof12_i),hb12a(tof12_i)) + 1
            tofmask(ch12b(tof12_i),hb12b(tof12_i)) =
     $           tofmask(ch12b(tof12_i),hb12b(tof12_i)) + 1

            tofmask(ch32a(tof32_i),hb32a(tof32_i)) =
     $           tofmask(ch32a(tof32_i),hb32a(tof32_i)) + 1
            tofmask(ch32b(tof32_i),hb32b(tof32_i)) =
     $           tofmask(ch32b(tof32_i),hb32b(tof32_i)) + 1

C-------

         ENDIF
      ENDIF

C     S21 - S31

      dist = ZTOF(3) - ZTOF(5)
      dl = 0.
      DO I=3,5
         dl = dl + TLOUT(i)
      ENDDO
      F = dl/dist

C      IF (tof21_i.GT.none_find.AND.tof31_i.GT.none_find) THEN
       IF ((tof21_i.GT.none_find).AND.(tof31_i.GT.none_find).AND.
     &    (xtofpos(2).NE.101.).AND.(ytofpos(3).NE.101.)) THEN
         t1 = tof21(1,tof21_i,itdc)
         t2 = tof21(2,tof21_i,itdc)
         t3 = tof31(1,tof31_i,itdc)
         t4 = tof31(2,tof31_i,itdc)
         IF ((t1.lt.4095).and.(t2.lt.4095).and.
     &        (t3.lt.4095).and.(t4.lt.4095)) THEN
            xhelp1 = tof21(1,tof21_i,itdc)+tof21(2,tof21_i,itdc)
            xhelp2 = tof31(1,tof31_i,itdc)+tof31(2,tof31_i,itdc)
            ds = xhelp1-xhelp2
            ihelp=(tof21_i-1)*3+tof31_i
            c1 = k_S21S31(1,ihelp)
            if (iz.gt.2) c1 = c1 + k1corrB1
            c2 = k_S21S31(2,ihelp)
            beta_a(5) = c2*F/(ds-c1)

C-------ToF Mask - S21 - S31

            tofmask(ch21a(tof21_i),hb21a(tof21_i)) =
     $           tofmask(ch21a(tof21_i),hb21a(tof21_i)) + 1
            tofmask(ch21b(tof21_i),hb21b(tof21_i)) =
     $           tofmask(ch21b(tof21_i),hb21b(tof21_i)) + 1

            tofmask(ch31a(tof31_i),hb31a(tof31_i)) =
     $           tofmask(ch31a(tof31_i),hb31a(tof31_i)) + 1
            tofmask(ch31b(tof31_i),hb31b(tof31_i)) =
     $           tofmask(ch31b(tof31_i),hb31b(tof31_i)) + 1

C-------

         ENDIF
      ENDIF

C     S21 - S32

      dist = ZTOF(3) - ZTOF(6)
      dl = 0.
      DO I=3,6
         dl = dl + TLOUT(i)
      ENDDO
      F = dl/dist

C      IF (tof21_i.GT.none_find.AND.tof32_i.GT.none_find) THEN
       IF ((tof21_i.GT.none_find).AND.(tof32_i.GT.none_find).AND.
     &    (xtofpos(2).NE.101.).AND.(xtofpos(3).NE.101.)) THEN
         t1 = tof21(1,tof21_i,itdc)
         t2 = tof21(2,tof21_i,itdc)
         t3 = tof32(1,tof32_i,itdc)
         t4 = tof32(2,tof32_i,itdc)
         IF ((t1.lt.4095).and.(t2.lt.4095).and.
     &        (t3.lt.4095).and.(t4.lt.4095)) THEN
            xhelp1 = tof21(1,tof21_i,itdc)+tof21(2,tof21_i,itdc)
            xhelp2 = tof32(1,tof32_i,itdc)+tof32(2,tof32_i,itdc)
            ds = xhelp1-xhelp2
            ihelp=(tof21_i-1)*3+tof32_i
            c1 = k_S21S32(1,ihelp)
            if (iz.gt.2) c1 = c1 + k1corrB1
            c2 = k_S21S32(2,ihelp)
            beta_a(6) = c2*F/(ds-c1)

C-------ToF Mask - S21 - S32

            tofmask(ch21a(tof21_i),hb21a(tof21_i)) =
     $           tofmask(ch21a(tof21_i),hb21a(tof21_i)) + 1
            tofmask(ch21b(tof21_i),hb21b(tof21_i)) =
     $           tofmask(ch21b(tof21_i),hb21b(tof21_i)) + 1

            tofmask(ch32a(tof32_i),hb32a(tof32_i)) =
     $           tofmask(ch32a(tof32_i),hb32a(tof32_i)) + 1
            tofmask(ch32b(tof32_i),hb32b(tof32_i)) =
     $           tofmask(ch32b(tof32_i),hb32b(tof32_i)) + 1

C-------

         ENDIF
      ENDIF

C     S22 - S31

      dist = ZTOF(4) - ZTOF(5)
      dl = 0.
      DO I=4,5
         dl = dl + TLOUT(i)
      ENDDO
      F = dl/dist
      
C WM workaround
      dl = dl - 0.06*F
      F = dl/dist

C      IF (tof22_i.GT.none_find.AND.tof31_i.GT.none_find) THEN
       IF ((tof22_i.GT.none_find).AND.(tof31_i.GT.none_find).AND.
     &    (ytofpos(2).NE.101.).AND.(ytofpos(3).NE.101.)) THEN
         t1 = tof22(1,tof22_i,itdc)
         t2 = tof22(2,tof22_i,itdc)
         t3 = tof31(1,tof31_i,itdc)
         t4 = tof31(2,tof31_i,itdc)
         IF ((t1.lt.4095).and.(t2.lt.4095).and.
     &        (t3.lt.4095).and.(t4.lt.4095)) THEN
            xhelp1 = tof22(1,tof22_i,itdc)+tof22(2,tof22_i,itdc)
            xhelp2 = tof31(1,tof31_i,itdc)+tof31(2,tof31_i,itdc)
            ds = xhelp1-xhelp2
            ihelp=(tof22_i-1)*3+tof31_i
            c1 = k_S22S31(1,ihelp)
            if (iz.gt.2) c1 = c1 + k1corrB1
            c2 = k_S22S31(2,ihelp)
            beta_a(7) = c2*F/(ds-c1)

C-------ToF Mask - S22 - S31

            tofmask(ch22a(tof22_i),hb22a(tof22_i)) =
     $           tofmask(ch22a(tof22_i),hb22a(tof22_i)) + 1
            tofmask(ch22b(tof22_i),hb22b(tof22_i)) =
     $           tofmask(ch22b(tof22_i),hb22b(tof22_i)) + 1

            tofmask(ch31a(tof31_i),hb31a(tof31_i)) =
     $           tofmask(ch31a(tof31_i),hb31a(tof31_i)) + 1
            tofmask(ch31b(tof31_i),hb31b(tof31_i)) =
     $           tofmask(ch31b(tof31_i),hb31b(tof31_i)) + 1

C-------

         ENDIF
      ENDIF

C     S22 - S32

      dist = ZTOF(4) - ZTOF(6)
      dl = 0.
      DO I=4,6
         dl = dl + TLOUT(i)
      ENDDO
      F = dl/dist
      
C WM workaround      
      dl = dl - 0.06*F
      F = dl/dist


C      IF (tof22_i.GT.none_find.AND.tof32_i.GT.none_find) THEN
       IF ((tof22_i.GT.none_find).AND.(tof32_i.GT.none_find).AND.
     &    (ytofpos(2).NE.101.).AND.(xtofpos(3).NE.101.)) THEN
         t1 = tof22(1,tof22_i,itdc)
         t2 = tof22(2,tof22_i,itdc)
         t3 = tof32(1,tof32_i,itdc)
         t4 = tof32(2,tof32_i,itdc)
         IF ((t1.lt.4095).and.(t2.lt.4095).and.
     &        (t3.lt.4095).and.(t4.lt.4095)) THEN
            xhelp1 = tof22(1,tof22_i,itdc)+tof22(2,tof22_i,itdc)
            xhelp2 = tof32(1,tof32_i,itdc)+tof32(2,tof32_i,itdc)
            ds = xhelp1-xhelp2
            ihelp=(tof22_i-1)*3+tof32_i
            c1 = k_S22S32(1,ihelp)
            if (iz.gt.2) c1 = c1 + k1corrB1
            c2 = k_S22S32(2,ihelp)
            beta_a(8) = c2*F/(ds-c1)

C-------ToF Mask - S22 - S32

            tofmask(ch22a(tof22_i),hb22a(tof22_i)) =
     $           tofmask(ch22a(tof22_i),hb22a(tof22_i)) + 1
            tofmask(ch22b(tof22_i),hb22b(tof22_i)) =
     $           tofmask(ch22b(tof22_i),hb22b(tof22_i)) + 1

            tofmask(ch32a(tof32_i),hb32a(tof32_i)) =
     $           tofmask(ch32a(tof32_i),hb32a(tof32_i)) + 1
            tofmask(ch32b(tof32_i),hb32b(tof32_i)) =
     $           tofmask(ch32b(tof32_i),hb32b(tof32_i)) + 1

C-------

         ENDIF
      ENDIF

C     S11 - S21

      dist = ZTOF(1) - ZTOF(3)
      dl = 0.
      DO I=1,3
         dl = dl + TLOUT(i)
      ENDDO
      F = dl/dist

C WM workaround      
      dl = dl - 0.442*F
      F = dl/dist

C      IF (tof11_i.GT.none_find.AND.tof21_i.GT.none_find) THEN
       IF ((tof11_i.GT.none_find).AND.(tof21_i.GT.none_find).AND.
     &    (ytofpos(1).NE.101.).AND.(xtofpos(2).NE.101.)) THEN
         t1 = tof11(1,tof11_i,itdc)
         t2 = tof11(2,tof11_i,itdc)
         t3 = tof21(1,tof21_i,itdc)
         t4 = tof21(2,tof21_i,itdc)
         IF ((t1.lt.4095).and.(t2.lt.4095).and.
     &        (t3.lt.4095).and.(t4.lt.4095)) THEN
            xhelp1 = tof11(1,tof11_i,itdc)+tof11(2,tof11_i,itdc)
            xhelp2 = tof21(1,tof21_i,itdc)+tof21(2,tof21_i,itdc)
            ds = xhelp1-xhelp2
            ihelp=(tof11_i-1)*2+tof21_i
            c1 = k_S11S21(1,ihelp)
            if (iz.gt.2) c1 = c1 + k1corrC1
            c2 = k_S11S21(2,ihelp)
            beta_a(9) = c2*F/(ds-c1)

C-------ToF Mask - S11 - S21

            tofmask(ch11a(tof11_i),hb11a(tof11_i)) =
     $           tofmask(ch11a(tof11_i),hb11a(tof11_i)) + 1
            tofmask(ch11b(tof11_i),hb11b(tof11_i)) =
     $           tofmask(ch11b(tof11_i),hb11b(tof11_i)) + 1

            tofmask(ch21a(tof21_i),hb21a(tof21_i)) =
     $           tofmask(ch21a(tof21_i),hb21a(tof21_i)) + 1
            tofmask(ch21b(tof21_i),hb21b(tof21_i)) =
     $           tofmask(ch21b(tof21_i),hb21b(tof21_i)) + 1

C-------

         ENDIF
      ENDIF

C     S11 - S22

      dist = ZTOF(1) - ZTOF(4)
      dl = 0.
      DO I=1,4
         dl = dl + TLOUT(i)
      ENDDO
      F = dl/dist

C      IF (tof11_i.GT.none_find.AND.tof22_i.GT.none_find) THEN
       IF ((tof11_i.GT.none_find).AND.(tof22_i.GT.none_find).AND.
     &    (ytofpos(1).NE.101.).AND.(ytofpos(2).NE.101.)) THEN
         t1 = tof11(1,tof11_i,itdc)
         t2 = tof11(2,tof11_i,itdc)
         t3 = tof22(1,tof22_i,itdc)
         t4 = tof22(2,tof22_i,itdc)
         IF ((t1.lt.4095).and.(t2.lt.4095).and.
     &        (t3.lt.4095).and.(t4.lt.4095)) THEN
            xhelp1 = tof11(1,tof11_i,itdc)+tof11(2,tof11_i,itdc)
            xhelp2 = tof22(1,tof22_i,itdc)+tof22(2,tof22_i,itdc)
            ds = xhelp1-xhelp2
            ihelp=(tof11_i-1)*2+tof22_i
            c1 = k_S11S22(1,ihelp)
            if (iz.gt.2) c1 = c1 + k1corrC1
            c2 = k_S11S22(2,ihelp)
            beta_a(10) = c2*F/(ds-c1)

C-------ToF Mask - S11 - S22

            tofmask(ch11a(tof11_i),hb11a(tof11_i)) =
     $           tofmask(ch11a(tof11_i),hb11a(tof11_i)) + 1
            tofmask(ch11b(tof11_i),hb11b(tof11_i)) =
     $           tofmask(ch11b(tof11_i),hb11b(tof11_i)) + 1

            tofmask(ch22a(tof22_i),hb22a(tof22_i)) =
     $           tofmask(ch22a(tof22_i),hb22a(tof22_i)) + 1
            tofmask(ch22b(tof22_i),hb22b(tof22_i)) =
     $           tofmask(ch22b(tof22_i),hb22b(tof22_i)) + 1

C-------

         ENDIF
      ENDIF

C     S12 - S21

      dist = ZTOF(2) - ZTOF(3)
      dl = 0.
      DO I=2,3
         dl = dl + TLOUT(i)
      ENDDO
      F = dl/dist

C  WM workaround
      dl = dl - 0.442*F
      F = dl/dist

C      IF (tof12_i.GT.none_find.AND.tof21_i.GT.none_find) THEN
       IF ((tof12_i.GT.none_find).AND.(tof21_i.GT.none_find).AND.
     &    (xtofpos(1).NE.101.).AND.(xtofpos(2).NE.101.)) THEN
         t1 = tof12(1,tof12_i,itdc)
         t2 = tof12(2,tof12_i,itdc)
         t3 = tof21(1,tof21_i,itdc)
         t4 = tof21(2,tof21_i,itdc)
         IF ((t1.lt.4095).and.(t2.lt.4095).and.
     &        (t3.lt.4095).and.(t4.lt.4095)) THEN
            xhelp1 = tof12(1,tof12_i,itdc)+tof12(2,tof12_i,itdc)
            xhelp2 = tof21(1,tof21_i,itdc)+tof21(2,tof21_i,itdc)
            ds = xhelp1-xhelp2
            ihelp=(tof12_i-1)*2+tof21_i
            c1 = k_S12S21(1,ihelp)
            if (iz.gt.2) c1 = c1 + k1corrC1
            c2 = k_S12S21(2,ihelp)
            beta_a(11) = c2*F/(ds-c1)

C-------ToF Mask - S12 - S21

            tofmask(ch12a(tof12_i),hb12a(tof12_i)) =
     $           tofmask(ch12a(tof12_i),hb12a(tof12_i)) + 1
            tofmask(ch12b(tof12_i),hb12b(tof12_i)) =
     $           tofmask(ch12b(tof12_i),hb12b(tof12_i)) + 1

            tofmask(ch21a(tof21_i),hb21a(tof21_i)) =
     $           tofmask(ch21a(tof21_i),hb21a(tof21_i)) + 1
            tofmask(ch21b(tof21_i),hb21b(tof21_i)) =
     $           tofmask(ch21b(tof21_i),hb21b(tof21_i)) + 1

C-------

         ENDIF
      ENDIF

C     S12 - S22

      dist = ZTOF(2) - ZTOF(4)
      dl = 0.
      DO I=2,4
         dl = dl + TLOUT(i)
      ENDDO
      F = dl/dist

C      IF (tof12_i.GT.none_find.AND.tof22_i.GT.none_find) THEN
       IF ((tof12_i.GT.none_find).AND.(tof22_i.GT.none_find).AND.
     &    (xtofpos(1).NE.101.).AND.(ytofpos(2).NE.101.)) THEN
         t1 = tof12(1,tof12_i,itdc)
         t2 = tof12(2,tof12_i,itdc)
         t3 = tof22(1,tof22_i,itdc)
         t4 = tof22(2,tof22_i,itdc)
         IF ((t1.lt.4095).and.(t2.lt.4095).and.
     &        (t3.lt.4095).and.(t4.lt.4095)) THEN
            xhelp1 = tof12(1,tof12_i,itdc)+tof12(2,tof12_i,itdc)
            xhelp2 = tof22(1,tof22_i,itdc)+tof22(2,tof22_i,itdc)
            ds = xhelp1-xhelp2
            ihelp=(tof12_i-1)*2+tof22_i
            c1 = k_S12S22(1,ihelp)
            if (iz.gt.2) c1 = c1 + k1corrC1
            c2 = k_S12S22(2,ihelp)
            beta_a(12) = c2*F/(ds-c1)

C-------ToF Mask - S12 - S22

            tofmask(ch12a(tof12_i),hb12a(tof12_i)) =
     $           tofmask(ch12a(tof12_i),hb12a(tof12_i)) + 1
            tofmask(ch12b(tof12_i),hb12b(tof12_i)) =
     $           tofmask(ch12b(tof12_i),hb12b(tof12_i)) + 1

            tofmask(ch22a(tof22_i),hb22a(tof22_i)) =
     $           tofmask(ch22a(tof22_i),hb22a(tof22_i)) + 1
            tofmask(ch22b(tof22_i),hb22b(tof22_i)) =
     $           tofmask(ch22b(tof22_i),hb22b(tof22_i)) + 1

C-------

         ENDIF
      ENDIF

C-------
C
C      icount=0
C      sw=0.
C      sxw=0.
C      beta_mean=100.
C
C      do i=1,12
C         if ((beta_a(i).gt.-1.5).and.(beta_a(i).lt.1.5)) then
C            icount= icount+1
C            if (i.le.4) w_i=1./(0.13**2.)
C            if ((i.ge.5).and.(i.le.8)) w_i=1./(0.16**2.)
C            if (i.ge.9) w_i=1./(0.25**2.) ! to be checked
C            sxw=sxw + beta_a(i)*w_i
C            sw =sw + w_i
C         endif
C      enddo
C
C      if (icount.gt.0) beta_mean=sxw/sw
C      beta_a(13) = beta_mean
C

C-------  New mean beta  calculation

       do i=1,12
         btemp(i) =  beta_a(i) 
       enddo

       beta_a(13)=newbeta(2,btemp,hitvec,10.,10.,20.)

C-------


c       IF (tof11_i.GT.none_find)
c     &  write(*,*) '11 ',tof11(1,tof11_i,itdc),tof11(2,tof11_i,itdc)
c       IF (tof12_i.GT.none_find)
c     &  write(*,*) '12 ',tof12(1,tof12_i,itdc),tof12(2,tof12_i,itdc)

c       IF (tof21_i.GT.none_find)
c     &  write(*,*) '21 ',tof21(1,tof21_i,itdc),tof21(2,tof21_i,itdc)
c       IF (tof22_i.GT.none_find)
c     &  write(*,*) '22 ',tof22(1,tof22_i,itdc),tof22(2,tof22_i,itdc)

c       IF (tof31_i.GT.none_find)
c     &  write(*,*) '31 ',tof31(1,tof31_i,itdc),tof31(2,tof31_i,itdc)
c       IF (tof32_i.GT.none_find)
c     &  write(*,*) '32 ',tof32(1,tof32_i,itdc),tof32(2,tof32_i,itdc)

c      write(*,*) xtofpos
c      write(*,*) ytofpos
C      write(*,*)'toftrk beta', beta_a
C      write(*,*) adcflagtof
C     write(*,*)'TOFTRK ',ievent,beta_a(1),beta_a(2),beta_a(3),beta_a(4)
c      write(*,*) 'toftrk'
c      write(*,*) xtofpos
c      write(*,*) ytofpos
c      write(*,*) xtr_tof
c      write(*,*) ytr_tof

c       write(*,*) '---------  end toftrk ----------'

      RETURN
      END


C------------------------------------------------------------------

