module OW_parms_areal

use sep
use OW_type

implicit none

logical,public, save				:: encode,adj
type (axis) ,public, save				:: kxm, h, w, x, xm, z, sh
type (dat_grid)	 ,public, save			:: rec, sou
type (mod_grid)	 ,public, save			:: img
type (velo) ,public, save				:: v
integer	 ,public, save				:: nref,iwb,tapertr,padtraces
integer	 ,public, save				:: isp,nwbott,node
real ,public, save					:: xmin,xmax,ymin,ymax,zmin,zmax,pi,dwbott
real, public, save					:: eps

contains

subroutine owparam(migrate)

logical						:: migrate
real						:: wbott,dd1,dd2,od1,od2,dv1,dv2,ov1,ov2,dzref,ozref
integer						:: iw,nt, i, m,nd1,nd2,nv1,nv2,nzref

if (migrate) then
call from_param("adj",adj,.true.)
if (adj) then
 call from_history("n1",rec%h%n);	call from_history("d1",rec%h%d);	call from_history("o1",rec%h%o)
 call from_history("n3",rec%w%n);	call from_history("d3",rec%w%d);	call from_history("o3",rec%w%o)
 call from_history("n4",rec%x%n);	call from_history("d4",rec%x%d);	call from_history("o4",rec%x%o)  

 call from_aux("sou","n1",sou%h%n);	call from_aux("sou","d1",sou%h%d);	call from_aux("sou","o1",sou%h%o)
 call from_aux("sou","n3",sou%w%n);	call from_aux("sou","d3",sou%w%d);	call from_aux("sou","o3",sou%w%o)
 call from_aux("sou","n4",sou%x%n);	call from_aux("sou","d4",sou%x%d);	call from_aux("sou","o4",sou%x%o)  

 call from_aux("vel","n1",v%z%n);	call from_aux("vel","d1",v%z%d);	call from_aux("vel","o1",v%z%o)
 call from_aux("vel","n2",v%xm%n);	call from_aux("vel","d2",v%xm%d);	call from_aux("vel","o2",v%xm%o)

 call from_aux("deltaS","n1",nd1);	call from_aux("deltaS","d1",dd1);	call from_aux("deltaS","o1",od1)
 call from_aux("deltaS","n2",nd2);	call from_aux("deltaS","d2",dd2);	call from_aux("deltaS","o2",od2)
else
 call from_history("n1",img%z%n);	call from_history("d1",img%z%d);	call from_history("o1",img%z%o)
 call from_history("n2",img%xm%n);	call from_history("d2",img%xm%d);	call from_history("o2",img%xm%o)  
 call from_history("n3",rec%w%n);	call from_history("d3",rec%w%d);	call from_history("o3",rec%w%o)
 call from_history("n4",rec%x%n);	call from_history("d4",rec%x%d);	call from_history("o4",rec%x%o)

 call from_aux("sou","n1",img%z%n);	call from_aux("sou","d1",img%z%d);	call from_aux("sou","o1",img%z%o)
 call from_aux("sou","n2",img%xm%n);	call from_aux("sou","d2",img%xm%d);	call from_aux("sou","o2",img%xm%o)  
 call from_aux("sou","n3",sou%w%n);	call from_aux("sou","d3",sou%w%d);	call from_aux("sou","o3",sou%w%o)
 call from_aux("sou","n4",sou%x%n);	call from_aux("sou","d4",sou%x%d);	call from_aux("sou","o4",sou%x%o)

 call from_aux("vel","n1",v%z%n);	call from_aux("vel","d1",v%z%d);	call from_aux("vel","o1",v%z%o)
 call from_aux("vel","n2",v%xm%n);	call from_aux("vel","d2",v%xm%d);	call from_aux("vel","o2",v%xm%o)

 call from_aux("deltaI","n1",img%z%n);	call from_aux("deltaI","d1",img%z%d);	call from_aux("deltaI","o1",img%z%o)
 call from_aux("deltaI","n2",img%sh%n);	call from_aux("deltaI","d2",img%sh%d);	call from_aux("deltaI","o2",img%sh%o)
 call from_aux("deltaI","n3",img%xm%n);	call from_aux("deltaI","d3",img%xm%d);	call from_aux("deltaI","o3",img%xm%o)  
end if
 call from_aux("vref","n1",nref)
 call from_aux("vref","n3",nzref);	call from_aux("vref","d3",dzref);	call from_aux("vref","o3",ozref)

 if (nzref .ne. v%z%n) then
    call seperr("# depth steps vref different from vel file")
 end if
 if (dzref .ne. v%z%d) then
    call seperr("depth interval vref different from vel file")
 end if
 if (ozref .ne. v%z%o) then
    call seperr("depth origin vref different from vel file")
 end if

 if (sou%w%n .ne. rec%w%n) then
    call seperr("# source freq. .ne. # receiver freq.")
 end if
 if (sou%x%n .ne. rec%x%n) then
    call seperr("# source planes. .ne. # receiver planes")
 end if
 if (sou%h%n .ne. rec%h%n) then
    call seperr("# source traces. .ne. # areal shot traces")
 end if
 if (sou%w%d .ne. rec%w%d) then
    call seperr("source freq.interval .ne. receiver freq.interval")
 end if
 if (sou%h%d .ne. rec%h%d) then
    call seperr("source trace interval .ne. areal shot trace interval")
 end if
 if (sou%x%d .ne. rec%x%d) then
    call seperr("source plane interval .ne. areal shot plane interval")
 end if
 if (sou%w%o .ne. rec%w%o) then
    call seperr("source freq.origin .ne. receiver freq.origin")
 end if
 if (sou%x%o .ne. rec%x%o) then
    call seperr("source plane origin .ne. areal shot plane origin")
 end if
 if (sou%w%o .ne. rec%w%o) then
    call seperr("source origin .ne. areal shot origin")
 end if

if (adj) then
   if (nd1.ne.v%z%n) then
      call seperr("deltaS # depth samples .ne. velocity interval ")
   end if

   if (nd2.ne.v%xm%n) then
      call seperr("deltaS # traces .ne. velocity interval ")
   end if

   if (dd1.ne.v%z%d) then
      call seperr("deltaS depth interval .ne. velocity interval ")
   end if
  
   if (dd2.ne.v%xm%d) then
      call seperr("deltaS trace interval .ne. velocity interval ")
   end if
  
   if (od1.ne.v%z%o) then
      call seperr("deltaS depth origin .ne. velocity interval ")
   end if
  
   if (od2.ne.v%xm%o) then
      call seperr("deltaS x origin .ne. velocity interval ")
   end if
else
end if
 !!  get parameters
  call from_param("nrefvel",nref) 
  call from_param("xmin",xmin) 
  call from_param("xmax",xmax)
  call from_param("ymin",ymin,0.)
  call from_param("ymax",ymax,0.)
  call from_param("zmin",zmin,0.)
  call from_param("zmax",zmax)
  call from_param("node",node,9999)


  if (xmin < v%xm%o ) call seperr("Min.x-image space < Min.x-velocity")
  if (xmax > v%xm%o+(v%xm%n-1)*v%xm%d ) call seperr("Max.x-image space > Max.x-velocity")
 
  if (zmin < v%xm%o ) call seperr("Min.z-image space < Min.z-velocity")
  if (zmax > v%z%o+(v%z%n-1)*v%z%d ) call seperr("Max.z-image space > Max.z-velocity")

  call from_param("wbott",wbott,0.)
  if(adj) then
     call from_param("nh",img%sh%n,1)
     call from_param("oh",img%sh%o,0.)
     call from_param("dh",img%sh%d,1.)
  end if
  call from_param("wbott",wbott,0.)
  call from_param("traces_in_taper",tapertr,0)
  call from_param("padtraces",padtraces,tapertr)
  call from_param("encode",encode,.false.)

!----------------------------------------------
!
! Migration output
!
if (adj) then
  img%z%n=floor((zmax-zmin)/v%z%d)+1
  img%z%o=zmin
  img%z%d=v%z%d
end if

if (adj) then
 call to_history("n1",img%z%n) ; call to_history("d1",img%z%d) ; call to_history("o1",img%z%o)
 call to_history("n2",img%sh%n) ; call to_history("d2",img%sh%d) ; call to_history("o2",img%sh%o)
 call to_history("n3",rec%h%n); call to_history("d3",rec%h%d); call to_history("o3",rec%h%o)  
 call to_history("n4",1); call to_history("d4",1); call to_history("o4",0)
 call to_history("n5",1); call to_history("d5",1); call to_history("o5",0)  
 call to_history("esize",4)
else
 call to_history("n1",img%z%n) ; call to_history("d1",img%z%d) ; call to_history("o1",img%z%o)
 call to_history("n2",1) ; call to_history("d2",1.) ; call to_history("o2",0.)
 call to_history("n3",img%xm%n); call to_history("d3",img%xm%d); call to_history("o3",img%xm%o)  
 call to_history("n4",1); call to_history("d4",1); call to_history("o4",0)
 call to_history("n5",1); call to_history("d5",1); call to_history("o5",0)  
 call to_history("esize",4)
end if

if(exist_file("ws")) then
 call to_history("n1",img%z%n,"ws") ; call to_history("d1",img%z%d,"ws") ; call to_history("o1",img%z%o,"ws")
if (adj) then
 call to_history("n2",rec%h%n,"ws"); call to_history("d2",rec%h%d,"ws"); call to_history("o2",rec%h%o,"ws")  
else
 call to_history("n2",img%xm%n,"ws"); call to_history("d2",img%xm%d,"ws"); call to_history("o2",img%xm%o,"ws")  
end if
 call to_history("n3",rec%w%n,"ws"); call to_history("d3",rec%w%d,"ws"); call to_history("o3",rec%w%o,"ws")
 call to_history("n4",rec%x%n,"ws"); call to_history("d4",rec%x%d,"ws"); call to_history("o4",rec%x%o,"ws")
 call to_history("n5",1,"ws"); call to_history("d5",1.,"ws"); call to_history("o5",1.,"ws")  
 call to_history("esize",8,"ws")
end if

if(exist_file("wr")) then
 call to_history("n1",img%z%n,"wr") ; call to_history("d1",img%z%d,"wr") ; call to_history("o1",img%z%o,"wr")
if (adj) then
 call to_history("n2",rec%h%n,"wr"); call to_history("d2",rec%h%d,"wr"); call to_history("o2",rec%h%o,"wr")  
else
 call to_history("n2",img%xm%n,"wr"); call to_history("d2",img%xm%d,"wr"); call to_history("o2",img%xm%o,"wr")  
end if
 call to_history("n3",rec%w%n,"wr"); call to_history("d3",rec%w%d,"wr"); call to_history("o3",rec%w%o,"wr")
 call to_history("n4",rec%x%n,"wr"); call to_history("d4",rec%x%d,"wr"); call to_history("o4",rec%x%o,"wr")
 call to_history("n5",1,"wr"); call to_history("d5",1.,"wr"); call to_history("o5",1.,"wr")  
 call to_history("esize",8,"wr")
end if


if (adj) then
  img%xm%n=rec%h%n+2*padtraces+2*tapertr
  img%xm%d=rec%h%d
  img%xm%o=rec%h%o-(padtraces+tapertr)*img%xm%d
else
  img%xm%n=img%xm%n+2*padtraces+2*tapertr
  img%xm%o=img%xm%o-(padtraces+tapertr)*img%xm%d
end if

  wbott=int(wbott/img%z%d)*img%z%d
  nwbott=int(wbott/img%z%d)+1
  dwbott=wbott/2.

pi=4*atan(1.)
isp=abs(floor(rec%h%o/rec%h%d)) + 1
!--------------------------------------------------
! Depth vector
!
allocate(img%z%all(img%z%n))
do iwb = 1,img%z%n
   img%z%all(iwb) = img%z%o + (iwb - 1) * img%z%d
end do

!--------------------------------------------------
! Frequency vector
! no zero frequency component
!
rec%w%o   = 2*pi * rec%w%o
rec%w%d   = 2*pi * rec%w%d
allocate(rec%w%all(rec%w%n))
do iwb = 1,rec%w%n
   rec%w%all(iwb) = rec%w%o + (iwb - 1) * rec%w%d
end do

!--------------------------------------------------
! Wavenumber vector
!
kxm%n = img%xm%n
allocate(kxm%all(kxm%n))
kxm%all=0.
if(kxm%n/=1) then
  m = kxm%n/2
  do i=1,m+1
     kxm%all(i)=i-1
  end do
  do i=m+2,kxm%n
     kxm%all(i)=i-kxm%n-1
  end do
else
  kxm%all=0. !! make sure k=0 if n=1
end if

kxm%all = kxm%all / (kxm%n*img%xm%d)      !! 1/length
kxm%all = 2 * PI * kxm%all     !! radians


!----------------------------------------------------
write(0,*) "   PARM   "
write(0,*) "=========="
write(0,*) "PARAMETERS"
write(0,*) "=========="
write(0,*) //
write(0,*) "no.ref.velocities   = ",nref
write(0,*) "min.x-coord         = ",xmin
write(0,*) "max.x-coord         = ",xmax
write(0,*) "min.y-coord         = ",ymin
write(0,*) "max.y-coord         = ",ymax
write(0,*) "min.depth           = ",zmin
write(0,*) "max.depth           = ",zmax
write(0,*) "traces in taper     = ",tapertr
write(0,*) "padded traces       = ",padtraces
write(0,*) //
write(0,*) "============"
write(0,*) "HEADER WORDS"
write(0,*) "============"
write(0,*) //
write(0,*) "SHOTS"
write(0,*) "====="
write(0,*) "number of recv  = ", rec%h%n
write(0,*) "recv spacing    = ", rec%h%d
write(0,*) "min.recv offset = ", rec%h%o
write(0,*) "number of freq  = ", rec%w%n
write(0,*) "freq spacing    = ", rec%w%d/(2*pi)
write(0,*) "minimum freq    = ", rec%w%o/(2*pi)
write(0,*) "number of shots = ", rec%x%n
write(0,*) "shot spacing    = ", rec%x%d
write(0,*) "min.shot-x      = ", rec%x%o
write(0,*) //
write(0,*) "SOURCE FUNCTION"
write(0,*) "==============="
write(0,*) "no.of traces    = ", sou%h%n
write(0,*) "trace spacing   = ", sou%h%d
write(0,*) "min.trace-x     = ", sou%h%o
write(0,*) "number of freq  = ", sou%w%n
write(0,*) "freq spacing    = ", sou%w%d
write(0,*) "minimum freq    = ", sou%w%o
write(0,*) "no.of s.fnct    = ", sou%x%n
write(0,*) "s.fnct spacing  = ", sou%x%d
write(0,*) "min.s.fnct-x    = ", sou%x%o
write(0,*) //
write(0,*) "VELOCITY MODEL"
write(0,*) "=============="
write(0,*) "no.of traces    = ", v%xm%n
write(0,*) "trace spacing   = ", v%xm%d
write(0,*) "min.trace-x     = ", v%xm%o
write(0,*) "no.depth samples= ", v%z%n
write(0,*) "depth interval  = ", v%z%d
write(0,*) "minimum depth   = ", v%z%o
write(0,*) //
write(0,*) "IMAGE SPACE"
write(0,*) "==========="
write(0,*) "no.of traces    = ", img%xm%n
write(0,*) "trace spacing   = ", img%xm%d
write(0,*) "min.trace-x     = ", img%xm%o
write(0,*) "no.depth samples= ", img%z%n
write(0,*) "depth interval  = ", img%z%d
write(0,*) "minimum depth   = ", img%z%o
write(0,*) //
write(0,*) "WAVENUMBER"
write(0,*) "=========="
write(0,*) "no.wavenumbers  = ", kxm%n
write(0,*) "waveno.spacing  = ", kxm%d
write(0,*) "waveno.origin   = ", kxm%o
write(0,*) //
!write(0,*) kxm%all
end if
end subroutine

end module
