program OWmig2d_areal

use sep 
use OW_pspi_areal
use OW_fftw_areal
use OW_type
use OW_grabvel_areal
use OW_exec_time
use OW_parms_areal
!use OW_datum
use OW_taper_areal

implicit none

real						:: start, finish
integer						:: ixs,ino,ni,node,i,st
real, allocatable				:: velint(:,:),taper(:)
real, allocatable				:: vel(:,:),R_image(:,:,:)
complex, allocatable				:: recv(:,:),source(:,:,:),tfld(:)
complex, allocatable				:: sou_dat(:,:,:)
complex,allocatable				:: G1(:,:,:)

integer, external 				:: omp_get_thread_num,omp_get_num_threads,omp_get_num_procs,omp_set_num_threads

 call sep_init()
 call sep_begin_prog()

node=omp_get_num_procs()
!st=omp_set_num_threads(floor(node/8.))

write(0,*) 'node = ',node

start=0.;finish=0.
call owexec_time(start)
!---------------------------------------
! IO
!
 call owparam(.true.)

!----------------------------------------------------
! Memory Allocation
!
allocate(source(sou%h%n,sou%w%n,sou%x%n))
allocate(velint(v%z%n,v%xm%n))
allocate(vel(img%z%n,img%xm%n))
allocate(G1(img%z%n,img%sh%n,rec%h%n)); 
allocate(taper(rec%h%n))

call owexec_time(finish)
write(0,'(8x,a,f10.2,a)')' Init.phase - Total elapsed time:', (finish - start),' (sec)'

call owtaper(tapertr,taper)

call sreed("vel",velint,v%z%n*v%xm%n*4)
call sreed("sou",source,sou%h%n*sou%w%n*sou%x%n*8)

!--------------------------------------------------
!
!
write(0,*) "DOWNWARD CONTINUES SOURCE AND RECEIVER WAVEFIELDS"
write(0,*) "================================================="

open(unit=2,file='lixo')

if (allocated(recv)) deallocate(recv);allocate(recv(rec%h%n,rec%w%n))

!!--------------------------------------------------
!! Setting FFTW plans
!!
!if (allocated(tfld)) deallocate(tfld);allocate(tfld(rec%h%n));tfld=0
!write(0,*) "pre-computation of tfld"
!call owfftw(tfld)
!write(0,*) "computed tfld"
!deallocate(tfld)


!--------------------------------------------------
! Migrating
!
!-----------------------------------------------------------------------------------------------      

    recv=0;G1=0

    start=0.;finish=0. 
    call owexec_time(start)

!--------------------------------------------------
! Reading data

    call sreed("in",recv,8*rec%h%n*rec%w%n)

    call owgrab_vel(velint,vel,ni)

if (exist_file("velout")) then
   call to_history("n1",img%z%n,"velout");call to_history("o1",v%z%o,"velout");call to_history("d1",v%z%d,"velout"); 
   call to_history("n2",img%xm%n,"velout");call to_history("o2",img%xm%o,"velout");call to_history("d2",img%xm%d,"velout");
   call srite("velout",vel,4*v%z%n*img%xm%n)
end if

!    if (wbott.ne.0.) then
!       allocate(sou_dat(rec%h%n,sou%w%n,sou%x%n));sou_dat=0.
!       sou_dat(isp,:,:)=source;
!       call owdatum(isp,wbott,vel,recv,sou_dat)
!       call owpspi_init(isp,ni,ixs,wbott,node)
!       call owpspi (recv,sou_dat,vsal,vel,G1)
!       deallocate(sou_dat)
!    else
       call owpspi_init(ni,ixs,wbott,node)
       call owpspi (taper,recv,source,vsal,vel,G1)
!    end if

    call owexec_time(finish)

allocate(R_image(img%z%n,img%sh%n,img%xm%n))
R_image=real(G1)
call srite("out",R_image,4*img%z%n*img%sh%n*img%xm%n)

deallocate(R_image,recv,source,velint,vel)
close(2)
call sep_end_prog()
end program OWmig2d_areal
