module OW_combine

use OW_type
use OW_parms_areal
use OW_exec_time

implicit none

contains

!------------------------------------------------------
subroutine owcombine(intmp,inset,nodes,rand1,igap)

character*80					:: hor
real,allocatable					:: pick1(:,:)
interger						:: par1(2)
integer						:: igap
integer,optional					:: nodes
real,optional					:: rand1(:,:,:)
integer						:: iw,iz,inode,ih
integer						:: jxs,ixs
real						:: start,finish
complex 						:: inset(:,:,:)
real						:: intmp(:,:,:)
real, allocatable					:: weight(:,:,:), temp(:,:)
integer, external 					:: omp_get_thread_num,omp_get_num_threads

if (encode.and.(rand.or.exist_file("InCode"))) then
   call owexec_time(start)
   !$OMP parallel do private (iw,iz,ixs,jxs,inode) schedule(dynamic,1) num_threads(nodes)
   do iw = 1,sou%w%n-subf
      do ixs=igap,sog%xm%n,ngap
         jxs=ixs+int((sog%xm%o+sog%sh%o-xmin)/sog%xm%d)+0.5
         do iz=nzmax,nzmin,-1
            if ((iz-1)*sog%z%d > 0.8) then
               inset(iz-nzmin+1,jxs:jxs+sog%sh%n-1,iw)=inset(iz-nzmin+1,jxs:jxs+sog%sh%n-1,iw)+ &
                                                 cmplx(intmp(iz,:,ixs),0)*exp(cmplx(0,eps*rand1(iw,ixs,1)))
            else
               inset(iz-nzmin+1,jxs:jxs+sog%sh%n-1,iw)=inset(iz-nzmin+1,jxs:jxs+sog%sh%n-1,iw)+ &
                                                 cmplx(intmp(iz,:,ixs),0)*exp(cmplx(0,eps*rand1(iw,ixs,2)))
            end if
         end do
      end do
   end do
   !$OMP end parallel do
   call owexec_time(finish)
   write(0,6)  "FINISHED PHASE-ENCODING ",int(finish-start),'(s)'
6  format(A24,i5,A4)
else if (encode.and.gold.or.encode.and.sgnpi) then
   if (exist_file("horsplit")) then
      open(51,file=hor,status='old')
      allocate(pick1(sog%xm%n,sog%sh%n));pick1=9999.
      read(51,*,end=1) par1
      do n=1,par1(2)
         read(51,*,end=1) pick1(par1(1)+n-1,:)
      end do
      1 close(51)
   end if
   call owexec_time(start)
   if (encz) then
      !$OMP parallel do private (ih,iw,&
                                   ixs,&
                                   jxs,&
                                  inode) schedule(dynamic,1) num_threads(nodes)
      do iw = 1,sou%w%n-subf
         do ixs=igap,sog%xm%n,ngap
            jxs=ixs+int((sog%xm%o+sog%sh%o-xmin)/sog%xm%d)+0.5
            do ih=1,sog%sh%n
            do iz=nzmax,nzmin,-1
!               if ((iz-1)*sog%z%d > 2000) then
               if (iz > pick1(jxs,ih)) then
                  inset(iz-nzmin+1,jxs:jxs+sog%sh%n-1,iw)=inset(iz-nzmin+1,jxs:jxs+sog%sh%n-1,iw)+ &
                                                    cmplx(intmp(iz,:,ixs),0)*exp(cmplx(0,eps*rand1(iw,1,1)))
!                                                    cmplx(intmp(iz,:,ixs),0)*exp(cmplx(0,eps*rand1(iw,ixs,1)))
               else
                  inset(iz-nzmin+1,jxs:jxs+sog%sh%n-1,iw)=inset(iz-nzmin+1,jxs:jxs+sog%sh%n-1,iw)+ &
                                                    cmplx(intmp(iz,:,ixs),0)*exp(cmplx(0,eps*rand1(iw,1,2)))
!                                                    cmplx(intmp(iz,:,ixs),0)*exp(cmplx(0,eps*rand1(iw,ixs,2)))
               end if
            end do
         end do
         end do
      end do
      !$OMP end parallel do
   else
      !$OMP parallel do private (   iw,&
                                   ixs,&
                                   jxs,&
                                  inode) schedule(dynamic,1) num_threads(nodes)
      do iw = 1,sou%w%n-subf
         do ixs=igap,sog%xm%n,ngap
            jxs=ixs+int((sog%xm%o+sog%sh%o-xmin)/sog%xm%d)+0.5
            inset(:,jxs:jxs+sog%sh%n-1,iw)=inset(:,jxs:jxs+sog%sh%n-1,iw)+ & 
                                  cmplx(intmp(nzmin:nzmax,:,ixs),0)*exp(cmplx(0,eps*rand1(iw,ixs,1)))
         end do
      end do
      !$OMP end parallel do

   end if
   call owexec_time(finish)
   write(0,6)  "FINISHED PHASE-ENCODING ",int(finish-start),'(s)'

else if (.not.encode) then
   do ixs=igap,sog%xm%n,ngap
      jxs=ixs+int((sog%xm%o+sog%sh%o-xmin)/sog%xm%d)+0.5
      inset(:,jxs:jxs+sog%sh%n-1,1)=inset(:,jxs:jxs+sog%sh%n-1,1)+cmplx(intmp(nzmin:nzmax,:,ixs),0)
   end do
end if
end subroutine 
end module

!if (exist_file("weight")) then
!   allocate(weight(sog%xm%n,sog%z%n,sou%w%n-subf),temp(sog%xm%n,sog%z%n))
!   call sreed("weight",temp,4*sog%xm%n*sog%z%n)
!   do iw=1,sou%w%n-subf
!      weight(:,:,iw)=temp      
!   end do
!   deallocate(temp)
!end if
! .
! .
! .
!               if (allocated(weight)) then
!                  inset(iz-nzmin+1,jxs:jxs+sog%sh%n-1,iw)=inset(iz-nzmin+1,jxs:jxs+sog%sh%n-1,iw)+ &
!                                                    cmplx(intmp(iz,:,ixs),0)*exp(cmplx(0,eps*(weight(ixs,iz,iw))))!!*rand1(iw,ixs,iz)
!               else
!                  inset(iz-nzmin+1,jxs:jxs+sog%sh%n-1,iw)=inset(iz-nzmin+1,jxs:jxs+sog%sh%n-1,iw)+ &
!                                                    cmplx(intmp(iz,:,ixs),0)*exp(cmplx(0,eps*rand1(iw,ixs,iz)))
!               end if
! .
! .
!    if (allocated(weight)) deallocate(weight)
