program OWspline
use sep 
use base_spline_mod1
use tipo_mod
use parametros_mod
implicit none

double precision, allocatable  		      :: vel(:),vel0(:)
double precision                              :: aa,bb
real, allocatable  		              :: in(:)
integer                                       :: stat

!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
 call sep_init()
 call sep_begin_prog()
 call parametro()
 allocate(in(4*v%xm%n*v%ym%n*v%z%n));in=0.
 if (axk%n.ne.99999 .and. ayk%n.ne.99999 .and. azk%n.ne.99999) then
    allocate(vel(4*axk%n*ayk%n*azk%n));vel=0.
 end if
 allocate(vel0(4*v%xm%n*v%ym%n*v%z%n));vel0=0.
 call sreed ("vel",in,4*v%xm%n*v%ym%n*v%z%n)
 vel0=in
 deallocate(in)
 if (chop) then
    aa=(1-chophi)*maxval(vel0)
    bb=(1-choplo)*minval(vel0)
    where(vel0>aa) vel0=aa
    where(vel0<bb) vel0=bb
 end if
 if (not(secgr))  then
    if (expgrad.ne.1) vel0=sign(abs(vel0)**expgrad,vel0)
 end if
 if (axk%n.ne.99999 .and. ayk%n.ne.99999 .and. azk%n.ne.99999) then
    stat = base_spln3d_init(v%xm,v%ym,v%z,axk,ayk,azk)
    vel=dble(0.)
    stat = base_spln3d_oper(.true. , .false., vel, vel0)
    vel0=dble(0.)
    stat = base_spln3d_oper(.false. ,.false., vel, vel0)
 end if
 if (expgrad.ne.1) then
    if (secgr) then
       vel0=sign(abs(vel0)**(expgrad*2),vel0)
       if (axk%n.ne.99999 .and. ayk%n.ne.99999 .and. azk%n.ne.99999) then
          stat = base_spln3d_init(v%xm,v%ym,v%z,axk,ayk,azk)
          vel=dble(0.)
          stat = base_spln3d_oper(.true. , .false., vel, vel0)
          vel0=dble(0.)
          stat = base_spln3d_oper(.false. ,.false., vel, vel0)
       end if
    end if
 end if
 vel0=(vel0/maxval(abs(vel0)))*minvel
 call srite("vspl",sngl(vel0),4*v%z%n*v%ym%n*v%xm%n)

 if (allocated(vel)) deallocate(vel)
 if (allocated(vel0)) deallocate(vel0)
 call sep_end_prog()
 end program 
