      REAL FUNCTION matra()
*********************************************************
*                                                       *
* This file was generated by HUWFUN.                    *
*                                                       *
*********************************************************
*
*     Ntuple Id:      2    
*     Ntuple Title:   Pamela data
*     Creation:       13/07/2004 11.39.53
*
*********************************************************
*
      INCLUDE 'common.txt'
*
      vector pede(2,22,96)
      vector matmip(2,22,96)
*
*--   Enter user code here
*      
      l = 0
      sthit = 0.
      do i=1,2
         do j=1,22
            pre = 0
            do k=1,96
               if (mod(k+15,16).eq.0) pre = pre + 1
               if (base(i,j,pre).eq.0.) then
                  inf = (pre - 1) * 16 + 1
                  sup = pre * 16 
                  valore = 0.
                  test = 0.
                  do m = inf, sup
                     valore = valore + dexyc(i,j,m) 
                     test = test + dexy(i,j,m)
                  enddo
c                  bas = valore / 16.
                  bas = base(i,j,pre)
c                  valore2 = dexyc(i,j,k) - pede(i,j,k)
                  if (test.ne.0.) then
                     valore2 = dexy(i,j,k) - pede(i,j,k)
                  else
                     valore2 = dexyc(i,j,k) - pede(i,j,k)
                  endif
                  IF (i.eq.1.and.pre.eq.6.and.j.eq.17) then
                     call hf1(900+j,float(k),valore2)
c                     print *,'i j k pre, pede',i,j,k,pre, pede(i,j,k)
cx                    print *,'base dexyc',base(i,j,pre), dexyc(i,j,k)
                  endif
               else
                  bas = base(i,j,pre)
                  IF (i.eq.1.and.pre.eq.6.and.j.eq.17) then
c                     valore2 = dexyc(i,j,k) - pede(i,j,k)
                     valore2 = dexy(i,j,k) - pede(i,j,k)
                     call hf1(900+j,float(k),valore2)
                     call hf1(999,float(k),base(i,j,pre))
                  endif
               endif
               if (dexyc(i,j,k).gt.0) then
c                  
                  sthit = sthit + 1.
                  calbase=(dexyc(i,j,k)-bas-pede(i,j,k))/matmip(i,j,k)
                  if (calbase.gt.0.7) then
                     if (i.eq.1) then
                        call hf2(300,float(k),float(j*2),calbase)
                     endif
                     if (i.eq.2) then
                        call hf2(400,float(k),float(j*2),calbase)
                     endif
                  endif
               endif
           enddo
         enddo
      enddo      
c
      call hf1(666,sthit,1.)
      matra = 1.
*
      END
