program s
#ifdef SEP_MPI
  use sep_par_mpi
#else
  use sep_par
#endif
  use patch_mod
  use srm_io_mod
  use verb_mod
  use semb_mod
  use off2ang_mod
  use  fft5d_mod
  implicit none
type(sep3d) :: in,out,gathers,wind
logical :: write_gathers,verb
character(len=128) :: intag,outtag,temp_tag
integer :: n(5),i,npatch,nown
integer, external :: sep_thread_num,sep3d_make_local_name
complex, pointer :: buf(:,:,:,:,:)
real, allocatable :: semb(:,:,:)
logical :: first_patch,first_gamma


#ifdef SEP_MPI
 call sep_par_start(intag=intag,outtag=outtag)
#endif
 call init_3d()
 call init_sep3d(intag,in,"INPUT")  !Initialize structure

 if(.not. init_verbose()) then
   call erexit("trouble initializing verbose")
 end if


 if (in%data_format(1:5).ne."FLOAT") then
   call erexit("only know how to handle float data")
 end if

  if(my_verb()) write(0,*)  "before initializing srm io"
 if(.not. init_srm_io(in,intag,out,outtag,gathers,"gathers",write_gathers)) &
   call seperr("trouble initializing output")


  if(my_verb()) write(0,*)  "before initializing patch"
 if(.not. init_patch_loop(in,wind,npatch)) &
   call seperr("trouble initializing patch loop")


  if(my_verb()) write(0,*)  "before initializing temp file"
if(0/=sep3d_make_local_name("tEMp",temp_tag)) then
    if(sep_thread_num()==0) write(0,*) "trouble make local temp file name"
    call seperr("")
end if


nown=0

do i=1,out%n(4)
  if(sep3d_own(out,i)) nown=nown+1
end do

if(.not. init_print_progress(0,npatch*nown)) &
  call seperr("trouble initializing print progress")

first_patch=.true.
  if(my_verb()) write(0,*)  "before looping patches"
do while(.not. done_set_patch(n,wind))
  allocate(buf(n(1),n(2),n(3),n(4),n(5)))
  allocate(semb(n(1),n(4),n(5)))

  if(my_verb()) write(0,*)  "before fft init"
  if(.not. fft_init(wind,buf,first_patch))&
    call seperr("trouble initializing fft")


  if(my_verb()) write(0,*)  "before read patch"
  if(.not. read_patch(intag,in,buf))  call seperr("trouble reading patch")


  if(my_verb()) write(0,*)  "before init srm"
  if(.not. init_srm(wind))&
   call seperr("trouble initializing srm")

  if(my_verb()) write(0,*)  "before init off2ang"
  if(.not. init_off2ang(wind))&
   call seperr("trouble initializing off2ang")

  if(my_verb()) write(0,*)  "before init semb"
  if(.not. init_semb(wind))&
    call seperr("trouble initializing semblance")

  if(my_verb()) write(0,*)  "before initial fft "
  if(.not. fft_forward(wind,buf))  &
    call seperr("trouble ffting data")


  call auxinout(temp_tag)
  if(.not. store_temp_patch(temp_tag,buf))  &
    call seperr("trouble storing patch")

  first_gamma=.true.
  do i=1,out%n(4) !depth,cmpx,cmpy,igam
    if(sep3d_own(out,i)) then

      if(.not.first_gamma) then
        if(.not. read_temp_patch(temp_tag,buf))  &
          call seperr("trouble reading patch")
      end if

      if(.not. srm_it(buf,(i-1)*out%d(4)+out%o(4))) &
        call seperr("trouble srming")

!    if(.not. sep3d_write_data("temp1.H",temp1,buf(:,:,:,:,1),&
!        nwind=(/gathers%n(1),gathers%n(2),gathers%n(3),gathers%n(4),gathers%n(5),1/),&
!        fwind=(/0,0,0,0,0,i-1/)))&
!           call seperr("trouble writing")
!
      if(.not. fft_backward(wind,buf))  &
        call seperr("trouble ffting data")
!    if(.not. sep3d_write_data("temp2.H",temp2,buf(:,:,:,:,1),&
!        nwind=(/gathers%n(1),gathers%n(2),gathers%n(3),gathers%n(4),gathers%n(5),1/),&
!        fwind=(/0,0,0,0,0,i-1/)))&
!           call seperr("trouble writing")
      if(.not. off2ang(buf))  &
        call seperr("trouble ffting data")

!    if(.not. sep3d_write_data("temp3.H",temp3,buf(:,:,:,:,1),&
!        nwind=(/gathers%n(1),gathers%n(2),gathers%n(3),gathers%n(4),gathers%n(5),1/),&
!        fwind=(/0,0,0,0,0,i-1/)))&
!           call seperr("trouble writing")
      if(.not. semb_it(buf,semb))  &
        call seperr("trouble ffting data")

      if(.not. write_patch(i,outtag,out,semb,"gathers",gathers,buf,write_gathers))  &
        call seperr("trouble ffting data")
      first_gamma=.false.
      call print_progress("finished srm pct")
    end if
  end do
  first_patch=.false.
  if(.not. srm_clean()) &
    call seperr("trouble cleaning srm")
  if(.not. off2ang_clean())&
    call seperr("trouble cleaning off2ang")
  if(.not. semb_clean())&
    call seperr("trouble cleaning semblance")
  deallocate(buf,semb)
end do
 call aux_unlink(temp_tag)
#ifdef SEP_MPI
 call sep_par_stop()
#endif
end program
