# 1 "<stdin>"
# 1 "<built-in>"
# 1 "<command line>"
# 1 "<stdin>"
!<
!interp_3d
!
!Usage:
!call interp_3d (arg1, arg2 arg3)
!
!Input Parameters:
!
!
!Output Parameters:
!
!
!
!Description:
!
!
!
!>
!-------------------------------------------------
!
!Author: Robert Clapp, ESMB 463, 7230253
!
!Date Created:Thu May 11 20:39:17 PDT 2000
!
!Purpose:
!
!
module put_3d_mod
  implicit none
  real,private, save :: o(3),d(3)
  integer,private, save :: n(3)
  contains
  subroutine init_coord(nn,oo,dd)
    integer :: nn(3)
    real :: oo(3),dd(3)
    n=nn
    o=oo
    d=dd
  end subroutine
  real function get_val_lin_pos(field,pos1,pos2,pos3)
    real :: field(:,:,:)
    real :: pos1,pos2,pos3
    real :: loc1,loc2,loc3
    integer :: iloc1,iloc2,iloc3
    loc1=(pos1-o(1))/d(1)+1.
    loc2=(pos2-o(2))/d(2)+1.
    loc3=(pos3-o(3))/d(3)+1.
    iloc1=loc1
    iloc2=loc2
    iloc3=loc3
    loc1=loc1-iloc1
    loc2=loc2-iloc2
    loc3=loc3-iloc3
    if (iloc1<1) then
      iloc1=1
      loc1=0.
    end if
    if (iloc1>=n(1)) then
      iloc1=n(1)-1
      loc1=1.
    end if
    if (iloc2<1) then
      iloc2=1
      loc2=0.
    end if
    if (iloc2>=n(2)) then
      iloc2=n(2)-1
      loc2=1.
    end if
    if (iloc3<1) then
      iloc3=1
      loc3=0.
    end if
    if (iloc3>=n(3)) then
      iloc3=n(3)-1
      loc3=1.
    end if
    if (n(3).eq.1) then
      get_val_lin_pos=( (1.-loc1)*(1.-loc2)*field(iloc1,iloc2,1)+(1.&
        &-loc1)*(loc2)*field(iloc1,iloc2+1,1)+(loc1)*(1.-loc2)*field&
        &(iloc1+1,iloc2,1)+(loc1)*(loc2)*field(iloc1+1,iloc2+1,1))
    else
      get_val_lin_pos=(1.-loc3)*( (1.-loc1)*(1.-loc2)*field(iloc1,iloc2&
        &,iloc3)+(1.-loc1)*(loc2)*field(iloc1,iloc2+1,iloc3)+(loc1)*&
        &(1.-loc2)*field(iloc1+1,iloc2,iloc3)+(loc1)*(loc2)*field(iloc1+1,iloc2+1,iloc3))
      get_val_lin_pos=get_val_lin_pos+(loc3)*( (1.-loc1)*(1.-loc2)*field&
        &(iloc1,iloc2,iloc3+1)+(1.-loc1)*(loc2)*field(iloc1,iloc2+1&
        &,iloc3+1)+(loc1)*(1.-loc2)*field(iloc1+1,iloc2,iloc3+1)+(loc1&
        &)*(loc2)*field(iloc1+1,iloc2+1,iloc3+1))
    end if
  end function
  subroutine get_interp_val(val,field,scale,iloc1,iloc2,iloc3,f1,f2,f3&
    &)
    real :: val,field(:,:,:),f1,f2,f3,scale
    integer :: iloc1,iloc2,iloc3
    if (n(3).eq.1) then
      val= val+(scale)*( (1.-f1)*(1.-f2)*field(iloc1,iloc2,1) +(f1)*&
        &(1.-f2)*field(iloc1+1,iloc2,1) +(1.-f1)*(f2)*field(iloc1&
        &,iloc2+1,1) +(f1)*(f2)*field(iloc1+1,iloc2+1,1))
    else
      val= val+(1.-f3)*(scale)*( (1.-f1)*(1.-f2)*field(iloc1,iloc2&
        &,iloc3) +(f1)*(1.-f2)*field(iloc1+1,iloc2,iloc3) +(1.-f1)*&
        &(f2)*field(iloc1,iloc2+1,iloc3) +(f1)*(f2)*field(iloc1+1,iloc2+1,iloc3))
      val= val+(f3)*(scale)*( (1.-f1)*(1.-f2)*field(iloc1,iloc2,iloc3&
        &+1)+(1.-f1)*(f2)*field(iloc1,iloc2+1,iloc3+1)+(f1)*(1.-f2)&
        &*field(iloc1+1,iloc2,iloc3+1)+(f1)*(f2)*field(iloc1+1,iloc2+1,iloc3+1))
    end if
  end subroutine
  subroutine put_interp_val(val,field,scale,iloc1,iloc2,iloc3,f1,f2,f3&
    &)
    real :: val,field(:,:,:),f1,f2,f3,scale
    integer :: iloc1,iloc2,iloc3
    if (n(3).eq.1) then
      field(iloc1,iloc2,1)= field(iloc1,iloc2,1)+scale*val*(1.-f1)*(1.&
        &-f2)
      field(iloc1+1,iloc2,1)= field(iloc1+1,iloc2,1)+scale*val*(f1)*&
        &(1.-f2)
      field(iloc1,iloc2+1,1)= field(iloc1,iloc2+1,1)+scale*val*(1.-f1&
        &)*(f2)
      field(iloc1+1,iloc2+1,1)= field(iloc1+1,iloc2+1,1)+scale*val*(f1&
        &)*(f2)
    else
      field(iloc1,iloc2,iloc3)= field(iloc1,iloc2,iloc3)+scale*val*(1.&
        &-f1)*(1.-f2)*(1.-f3)
      field(iloc1+1,iloc2,iloc3)= field(iloc1+1,iloc2,iloc3)+scale*val&
        &*(f1)*(1.-f2)*(1.-f3)
      field(iloc1,iloc2+1,iloc3)= field(iloc1,iloc2+1,iloc3)+scale*val&
        &*(1.-f1)*(f2)*(1.-f3)
      field(iloc1+1,iloc2+1,iloc3)= field(iloc1+1,iloc2+1,iloc3)+scale&
        &*val*(f1)*(f2)*(1.-f3)
      if (f3 > 0.01) then
        field(iloc1,iloc2,iloc3+1)= field(iloc1,iloc2,iloc3+1)+scale&
          &*val*(1.-f1)*(1.-f2)*(f3)
        field(iloc1+1,iloc2,iloc3+1)= field(iloc1+1,iloc2,iloc3+1)&
          &+scale*val*(f1)*(1.-f2)*(f3)
        field(iloc1,iloc2+1,iloc3+1)= field(iloc1,iloc2+1,iloc3+1)&
          &+scale*val*(1.-f1)*(f2)*(f3)
        field(iloc1+1,iloc2+1,iloc3+1)= field(iloc1+1,iloc2+1,iloc3&
          &+1)+scale*val*(f1)*(f2)*(f3)
      end if
    end if
  end subroutine
end module
