program OWmig2d_par_DD

use sep 
use OW_pspi_dd
use OW_fftw
use OW_type
use OW_grabvel
use OW_exec_time
use OW_parms
use OW_image
use OW_sumimg


implicit none

real						:: start, finish
integer						:: ixs,ino,ni,isp,inode,node,i
real, allocatable				:: velint(:,:)
real, allocatable				:: vel(:,:)
complex, allocatable				:: G1(:,:,:)
complex, allocatable				:: source(:,:,:),image(:,:,:)

integer, external 				:: omp_get_thread_num,omp_get_num_threads,omp_get_num_procs

 call sep_init()
 call sep_begin_prog()

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

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

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

!--------------------------------------------------
! Setting FFTW plans
!
!call owfftw(tfld)

!--------------------------------------------------
! Reading source fctn and velocity
 call sreed("sou",source,sou%h%n*sou%w%n*sou%x%n*8)
 call sreed("vel",velint,v%z%n*v%xm%n*4)

 source(1,:,1)=source(1,:,1)*cmplx(0,rec%w%all/(2*pi))

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

!--------------------------------------------------
!
!
write(0,*) "DOWNWARD CONTINUES SOURCES AND RECEIVERS"
write(0,*) "========================================"

!--------------------------------------------------
! shot position in the acq.geom.
isp = abs(floor(rec%h%o/rec%h%d)) + 1

open(unit=2,file='lixo')
!--------------------------------------------------
! Migrating
!
!-----------------------------------------------------------------------------------------------      
shot:do ixs=1,rec%x%n

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

! $OMP parallel default(shared) private(inode,ixs,vel,ni,GG)   

    inode=omp_get_thread_num()+1;

!    GG(:,:,:,inode+1)=0.

    call owgrab_vel(ixs,isp,velint,vel,ni)
    call owpspi_init(isp,ni,inode,node,ixs,500.)
    call owpspi (source,vel,image,G1)

! $OMP end parallel

!    call owsum_img_init(ni,node)
!    call owsum_img(GG(:,:,:,inode),image)

!    image(:,:,ni:ni+rec%h%n-1)=image(:,:,ni:ni+rec%h%n-1)+GG(:,:,:,inode+1);

    call owexec_time(finish)
    write(0,*) " DONE SHOT # ",ixs," OUT OF ",rec%x%n," in ", (finish - start)/60.," (min)"

end do shot

if(exist_file("G1")) then
 rec%w%d=rec%w%d/(2*pi); rec%w%o=rec%w%o/(2*pi);
 call to_history("n1",img%z%n,"G1") ; call to_history("d1",img%z%d,"G1") ; call to_history("o1",img%z%o,"G1")
 call to_history("n2",rec%h%n,"G1"); call to_history("d2",rec%h%d,"G1"); call to_history("o2",rec%h%o,"G1")
 call to_history("n3",rec%w%n,"G1"); call to_history("d3",rec%w%d,"G1"); call to_history("o3",rec%w%o,"G1")
 call to_history("esize",8,"G1")
 call srite("G1",G1,8*img%z%n*rec%h%n*rec%w%n)
end if

call srite("out",real(image),4*img%z%n*img%sh%n*img%xm%n)
close(2)
call sep_end_prog()
end program OWmig2d_par_DD
