!!.. read nzs depth step velocity or slowness
module wem_velref_mod
  implicit none
  integer, private :: nvrf
contains

!=====================================================================
!
  integer function wem_velref_init(nvrf_in) result(stat)
    integer :: nvrf_in
    nvrf  = nvrf_in
    stat  = 0
  end function
!=====================================================================
!
  integer function wem_velref_oper(vmod_slice, vref_slice) result(stat)
    real, dimension(:,:) :: vmod_slice
    real, dimension(:)   :: vref_slice
    stat = 1
    call wem_get_vrf_slice(vmod_slice, vref_slice)
    stat = 0
  end function
!======================================================================
  integer function wem_velref_clos() result(stat)
    !lloyd = .false.
    nvrf  = 0
    stat  = 0
  end function
!======================================================================
!.. select the reference velocity for the current depth step
!.. Lloyd's method will be used
!----------------------------------------------------------------------
  subroutine wem_get_vrf_slice(vmod_slice, vref_slice)
    real, dimension(:,:) :: vmod_slice
    real, dimension(:)   :: vref_slice
    integer              :: iref
    real                 :: vmax, vmin, dvel

    if(nvrf > size(vref_slice)) call seperr("ERROR: nref larger than the size of vref_slice(:)")

    !if (lloyd) then
    !  write(0,*)"Currently not available..."
    !  STOP
    !else
      vref_slice=0.0
      if(nvrf==1)then
         vref_slice(1)=sum(vmod_slice)/size(vmod_slice)
      else
         vmax=maxval(vmod_slice)
         vmin=minval(vmod_slice)
         dvel=(vmax-vmin)/(nvrf-1)
         do iref=1, nvrf
            vref_slice(iref)=vmin+(iref-1)*dvel
         end do
      end if
    !end if
  end subroutine

end module
