program OW_encode1

use sep 
use OW_xwell1
use OW_fftw
use OW_type
use OW_exec_time
use OW_parms
use OW_rndseed
use OW_hadam1
use OW_mseq
use wei_myfft_mod

implicit none

real						:: start, finish,a(2),b
integer						:: ixs,isp,i,stat,igap,n1,rxs,sxs,n_enc_sht,jgap,jxs,ireal,suba
integer, allocatable				:: CHD(:)
real, allocatable					:: velint(:,:),vel(:,:),codeX(:)
real, allocatable					:: rand1(:),code(:,:),taper(:),phase(:,:),hold(:,:),gol(:,:)
complex, allocatable				:: source(:,:,:),soutemp(:),tfld(:),wvfld(:,:),wvfld1(:,:,:)

 call sep_init()
 call sep_begin_prog()

!---------------------------------------
! IO
!
 call owparam()
 allocate(taper(sou%h%n))
 call owtaper(int(sou%h%n/10.),taper)
!----------------------------------------------------
! Memory Allocation
!
allocate(source(sou%h%n,sou%w%n,sou%x%n))
allocate(vel(img%z%n,kxm%n),velint(v%z%n,v%xm%n))
!--------------------------------------------------
! Reading source fctn and velocity
 call sreed("in",source,sou%w%n*sou%x%n*sou%h%n*8)
 call sreed("vel",velint,v%z%n*v%xm%n*4)
n1 = floor((sou%x%o-v%xm%o)/v%xm%d)+1
vel=velint(1:img%z%n,n1:n1+kxm%n-1);

deallocate(velint)
!--------------------------------------------------
!
!
if (random) then
   open(unit=4,file='random.txt')
else if (mseq) then
   open(unit=4,file='mseq.txt')
else if (gold) then
   open(unit=4,file='gold.txt')
end if
!--------------------------------------------------
! Encoding shots
!
!
start=0.;finish=0. 
call owexec_time(start)
if (not(noencode)) then

if (random) then
   allocate(code(nreal*sou%x%n,sou%w%n-subf))
   allocate(phase(nreal*sou%x%n,sou%w%n-subf))
   do ixs=1,nreal*sou%x%n
      allocate(rand1(sou%w%n-subf))
      call init_random_seed(ixs,sou%w%n-subf,rand1)
      code(ixs,:)=rand1
      deallocate(rand1)
   end do
   phase=cos(pi*code)
   suba=subf

else if(mseq) then
   allocate(code(nreal*sou%x%n,sou%w%n-subf))
   allocate(phase(nreal*sou%x%n,sou%w%n-subf))
   allocate(tfld(sou%w%n-subf))
!   call owfftw(tfld)
   do ixs=1,nreal*sou%x%n
      allocate(rand1(sou%w%n-subf));rand1=0
      call owmseq(2,8,rand1,2*(sou%w%n),mod(16+ixs,16))
      rand1=2*rand1-1
      tfld=0.;tfld=cmplx(rand1,0)
      phase(ixs,:)=rand1
!      call sfftw_execute(planF)
      call FFT1DT(tfld, 1)
      code(ixs,:)=atan(imag(tfld)*real(tfld)/(real(tfld)*real(tfld)+tiny(1.)))
      deallocate(rand1)
   end do
   deallocate(tfld)
   suba=subf

!else if(golay) then
!   allocate(code(2*sou%w%n,sou%w%n))
!   allocate(CHD(sou%w%n,sou%w%n))
!   call owhadamard(CHD)
!   do i=1,n_enc_sht
!      code(i,:)=CHD((i-1)*int(sou%w%n/n_enc_sht*0.5)+1,:)
!      code(i+n_enc_sht,:)=CHD((i-1)*int(sou%w%n/n_enc_sht*0.5)+int(sou%w%n/2.)+1,:)
!   end do
!   deallocate(CHD)

else if(gold) then
   allocate(code(nreal*sou%x%n,sou%w%n))
   allocate(hold(2,sou%w%n-subf))
   allocate(gol(sou%w%n,sou%w%n))
   allocate(phase(nreal*sou%x%n,sou%w%n))
!   allocate(tfld(sou%w%n))
!   call owfftw(tfld)
!!!!!!!!   do ixs=1,nreal*sou%x%n
   do ixs=1,2
      allocate(rand1(sou%w%n-subf));rand1=0
      call owmseq(2,8,rand1,2*(sou%w%n),mod(16+ixs,16))
      hold(ixs,:)=rand1
      deallocate(rand1)
   end do
   hold=hold*2-1
   do i=1,sou%w%n-subf
      gol(i,1:sou%w%n-subf)=hold(1,:)*cshift(hold(2,:),i)
   end do
   allocate(tfld(sou%w%n))
   do i=1,nreal*sou%x%n
      if (gol(i,1).eq.1) then
         gol(i,sou%w%n)=-1
      else
         gol(i,sou%w%n)=1
      end if
      phase(i,:)=gol(i,:)
      tfld=cmplx(gol(i,:),0)
      call FFT1DT(tfld, 1)
      code(i,:)=atan(imag(tfld)*real(tfld)/(real(tfld)*real(tfld)+tiny(1.)))
   end do
   deallocate(tfld,hold,gol)   
   suba=0
   

!!!!!!!!   allocate(CHD(sou%w%n-subf))
!!!!!!!!   do ixs=1,nreal*sou%x%n
!!!!!!!!!      if (ixs.le.sou%x%n) then
!!!!!!!!!         CHD=ieor(int(code(ixs,:)),int(code(ixs+(nreal-1)*sou%x%n,:)))
!!!!!!!!!      else
!!!!!!!!!         CHD=ieor(int(code(ixs,:)),int(code(ixs-1,:)))
!!!!!!!!!      end if
!!!!!!!!      if (ixs==1) then
!!!!!!!!         where(int(code(ixs,:)).eq.0.and.int(code(nreal*sou%x%n,:)).eq.0.or. &
!!!!!!!!               int(code(ixs,:)).eq.1.and.int(code(nreal*sou%x%n,:)).eq.1)
!!!!!!!!               CHD=1
!!!!!!!!         elsewhere
!!!!!!!!               CHD=0
!!!!!!!!         end where
!!!!!!!!      else
!!!!!!!!         where(int(code(ixs,:)).eq.0.and.int(code(ixs-1,:)).eq.0.or. &
!!!!!!!!               int(code(ixs,:)).eq.1.and.int(code(ixs-1,:)).eq.1)
!!!!!!!!               CHD=1
!!!!!!!!         elsewhere
!!!!!!!!               CHD=0
!!!!!!!!         end where
!!!!!!!!      end if

!!!!!!!!write(0,*) CHD

!!!!!!!!      tfld=cmplx(float(CHD),0)
!!!!!!!!      call FFT1DT(tfld, 1)
!!!!!!!!!      call sfftw_execute(planF)
!!!!!!!!      code(ixs,:)=atan(imag(tfld)*real(tfld)/(real(tfld)*real(tfld)+tiny(1.)))
!!!!!!!!   end do   
!!!!!!!!   deallocate(tfld,CHD)
else if(xiao) then
   allocate(codeX(nreal*sou%x%n))
   call init_random_seed(100,size(a),a)
   codeX(1)=a(1)
   do i=2,nreal*sou%x%n
      codeX(i)=atan(-(sum(cos(codex(1:i))/sin(codex(1:i)))))
   end do
   suba=0
end if
else
   allocate(code(nreal*sou%x%n,sou%w%n-subf))
   code=0.
   suba=subf

end if
call owexec_time(finish)
write(0,*) "Finished generating codes in ", finish-start," (s) ", suba

if(suba==0.and.subf.ne.2) then
   write(0,*) "vou escrever 0"
   do i=1,nreal*sou%x%n
        write(4,1111) phase(i,:) 
   1111 format(256f10.6)
   end do
else if(suba==1) then
   write(0,*) "vou escrever 1"
   do i=1,nreal*sou%x%n
        write(4,1112) phase(i,:)
   1112 format(255f10.6)
   end do
end if
if (outcode) then

write(0,*) nreal*sou%x%n,sou%x%d,sou%x%o

 call to_history("n1",nreal*sou%x%n,"code"); call to_history("d1",sou%x%d,"code"); call to_history("o1",sou%x%o,"code")
 call to_history("n2",1,"code"); call to_history("d2",0.,"code"); call to_history("o2",0.,"code")
 call to_history("n3",sou%w%n-suba,"code") ; call to_history("d3",sou%w%d/(2*pi),"code") ; call to_history("o3",sou%w%o/(2*pi),"code")
 call to_history("n4",1,"code"); call to_history("d4",1.,"code"); call to_history("o4",0.,"code")
 call to_history("esize",4,"code")
 call srite("code",code,4*nreal*sou%x%n*(sou%w%n-suba))
end if

allocate(soutemp(sou%w%n-suba))

do ixs=1,nreal*sou%x%n

   ireal=int((ixs-1)/sou%x%n*1.)+1
   jxs=ixs-int((ixs-1)/sou%x%n*1.)*sou%x%n
   if (xiao) then
      soutemp=source(401,1:sou%w%n-suba,1)*exp(cmplx(0,2*pi*codeX(ixs)))
   else
      soutemp=source(401,1:sou%w%n-suba,1)*exp(cmplx(0,2*pi*code(ixs,:)))
   end if
   call owxwell(ireal,jxs,soutemp,vel,taper)

end do

if (exist_file("decoded")) then
allocate(wvfld(sou%h%n,nreal*(sou%w%n-suba)),wvfld1(sou%h%n,sou%w%n-suba,sou%x%n))
wvfld=0.; wvfld1=0.
call sreed("wvfld",wvfld,8*sou%h%n*sou%w%n*nreal)
do ixs=1,nreal*sou%x%n
   jxs=ixs-int((ixs-1)/sou%x%n*1.)*sou%x%n
   do i=1,sou%h%n
      if (xiao) then
         wvfld1(i,:,jxs)=wvfld1(i,:,jxs)+wvfld(i,1:sou%w%n-suba)*exp(cmplx(0,-2*pi*codeX(ixs)))
      else
         wvfld1(i,:,jxs)=wvfld1(i,:,jxs)+wvfld(i,1:sou%w%n-suba)*exp(cmplx(0,-2*pi*code(ixs,:)))
      end if
   end do
end do
call srite("decoded",wvfld1,8*sou%x%n*sou%h%n*(sou%w%n-suba))
deallocate(wvfld,wvfld1)
end if

close(4)

call sep_end_prog()
end program OW_encode1
