# 1 "<stdin>"
# 1 "<built-in>"
# 1 "<command line>"
# 1 "<stdin>"
!<;
!dip_field
!
!Usage:
!dip_field.x <in.H >out.H
!
!Input Parameters:
!
!
!Output Parameters:
!
!
!
!Description:
!
!
!
!>
!-------------------------------------------------
!
!Author: Robert Clapp, ESMB 463, 7230253
!
! Wed Jan 8 14:52:30 PST 2003
!
!Purpose:
!
! Smooth a 3-D field
!
!
module smooth_it_mod
  implicit none
  contains
  logical function smooth_dips(ns,ds,sxz,syz,coher,nsmooth_dips)
    integer :: ns(3),nsmooth_dips(3)
    real :: ds(3)
    real :: sxz(:,:,:),syz(:,:,:),coher(:,:,:)
    real, allocatable :: pxz(:,:,:),pyz(:,:,:)
    real :: denom1,denom2,dy,dx,dz
    integer :: i1,i2,i3,b1,b2,b3,e1,e2,e3,ia,ib,ic
    allocate(pyz(ns(1),ns(2),ns(3)))
    allocate(pxz(ns(1),ns(2),ns(3)))
!SMOOTH THE DIPS
    pxz=0
    if (ns(3).ne.1) then
      pyz=0.
    end if
    do i3=1,ns(3)
      b3=max(1,i3-nsmooth_dips(3))
      e3=min(ns(3),i3+nsmooth_dips(3))
      do i2=1,ns(2)
        b2=max(1,i2-nsmooth_dips(2))
        e2=min(ns(2),i2+nsmooth_dips(2))
        do i1=1,ns(1)
          b1=max(1,i1-nsmooth_dips(1))
          e1=min(ns(1),i1+nsmooth_dips(1))
          denom1=0.
          denom2=0.
          do ic=b3,e3
            dy=(ic-i3+.2)*ds(3)
            do ib=b2,e2
              dx=(ib-i2+.2)*ds(2)
              do ia=b1,e1
                dz=(ia-i1+.2)*ds(1)
                pxz(i1,i2,i3)= pxz(i1,i2,i3)+sxz(ia,ib,ic)
                !*coher(ia,ib,ic)**3/(dx*dx+dz*dz+dy*dy)
                denom1= denom1+1
                !coher(ia,ib,ic)**3/(dx*dx+dz*dz+dy*dy)
                if (ns(3).ne.1) then
                  pyz(i1,i2,i3)= pyz(i1,i2,i3)+syz(ia,ib,ic)
                  !*coher(ia,ib,ic)**3/(dx*dx+dz*dz+dy*dy)
                  denom2= denom2+1
                  !coher(ia,ib,ic)**3/(dx*dx+dz*dz+dy*dy)
                end if
              end do
            end do
          end do
          if (denom1.ne.0.) then
            pxz(i1,i2,i3)=pxz(i1,i2,i3)/denom1
          end if
          if (ns(3).ne.1 .and. denom2.ne.0.) then
            pyz(i1,i2,i3)=pyz(i1,i2,i3)/denom2
          end if
        end do
      end do
    end do
    sxz=pxz
    if (ns(3).ne.1) then
      syz=pyz
    end if
    deallocate(pxz,pyz)
    smooth_dips=.true.
  end function
end module
