program OW_editvel_pikH

use sep

implicit none

logical			:: vmut,vper,vvel,vchg,neg
real			:: o1,o2,d1,d2,v,x
integer			:: i1,i2,n1,n2,shift
integer, allocatable          :: i(:)
real,allocatable		:: velin(:,:)

call sep_init()
call sep_begin_prog()

call from_history("n1",n2); call from_history("d1",d2); call from_history("o1",o2)
call from_history("n3",n1); call from_history("d3",d1); call from_history("o3",o1)

call to_history("n1",n2); call to_history("d1",d2); call to_history("o1",o2)
call to_history("n2",1); call to_history("d2",1.); call to_history("o2",0.)
call to_history("n3",1); call to_history("d3",1.); call to_history("o3",0.)
allocate(velin(n2,n1));velin=0;
allocate(i(n2));i=0
call sreed("in",velin,4*n1*n2)
call from_param("v",v)
call from_param("shift",shift,0)
if (shift.ne.0) call from_param("neg",neg,.false.)

   outer:do i2=1,n2
      do i1=1,n1
         if (velin(i2,i1)>=v) then
            if (shift.ne.0) then
               if (neg) then
                  i(i2)=i1-shift
               else
                  i(i2)=shift-i1+1
               end if
            else
               i(i2)=i1
            end if
         end if
      end do
   end do outer

call srite("out",float(i-1)*d1,4*n2)
deallocate(velin,i)
call sep_end_prog()
end program
