program run_it
  use io_nfilt
  use sep
  integer :: n(3),nf
  real :: o(3),oang,dang
  real :: d(3),slp
  character(len=128) :: label
  type(nfilter) :: filt
  real,allocatable :: dips(:),bank(:,:,:)
  integer :: nb,nang,iang,i,ncoef,j,ifilt
  integer,allocatable :: pch(:),n_s(:)
  call sep_init()
  call sep_get_data_axis_par("bank",1,ncoef,o(1),d(1),label)
  call sep_get_data_axis_par("bank",2,nang,oang,dang,label)
  call sep_get_data_axis_par("bank",3,nf,o(1),d(1),label)

  call sep_get_data_axis_par("in",1,n(1),o(1),d(1),label)
  call sep_get_data_axis_par("in",2,n(2),o(2),d(2),label)
  call sep_get_data_axis_par("in",3,n(3),o(3),d(3),label)

  allocate(pch(product(n)),dips(product(n)))
  allocate(bank(ncoef,nang,nf),n_s(nang))
  call sreed("in",dips,size(dips)*4)
  call sreed("bank",bank,size(bank)*4)
  call from_param("filt",ifilt,1)
  n_s=ncoef

  do i=1,product(n)
     slp=-45./atan(1.)*atan(dips(i))
     iang=max(1,min(nang,nint((slp-oang)/dang)))
     pch(i)=iang
   !  write(0,*) iang,i
  end do
  call nallocate(filt,n_s,pch)
  do i=1,nang
    do j=1,ncoef/2
      filt%hlx(i)%lag(j)=j
      filt%hlx(i)%lag(j+ncoef/2)=n(1)-ncoef/2+j+1
      filt%hlx(i)%flt(j)=bank(j,i,ifilt)
      filt%hlx(i)%flt(j+ncoef/2)=bank(j+ncoef/2,i,ifilt)
    end do
   ! write(0,*) i,filt%hlx(i)%lag
  end do
 ! call nallocate(filt,n_s,pch)
  
  if(.not. rite_nfilt("out",filt,n,o,d)) then
    call seperr("trouble writing out filter")
  end if
  
  
  
  
end program