program main
  use sep
  use off2ang2d_sd_transp_mod
  use solver_tiny_mod
  use cgstep_mod
  implicit none


  real, dimension(:), pointer :: modl(:), data(:)
  integer                     :: nz, nhx, nax, nmx, imx, stat, style
  real                        :: z0, hx0, mx0 
  real                        :: dz, dhx, dax, dmx
  real                        :: axmin, axmax 
  real, parameter             :: PI = 3.14159265359
  logical, parameter          :: T = .true. , F = .false.
  call initpar()

  call from_history("n1", nz ); call from_history("d1", dz ); call from_history("o1", z0 )
  call from_history("n2", nhx); call from_history("d2", dhx); call from_history("o2", hx0)  
  call from_history("n3", nmx); call from_history("d3", dmx); call from_history("o3", mx0)

  call from_param("amin", axmin, -60.)
  call from_param("amax", axmax, 60.)
  call from_param("style", style, 0)
  nax = nhx
  if (nax == 1) then
      dax = 0
  else
      dax = (axmax - axmin) / (nax - 1)
  end if

  allocate(modl(nz*nhx), data(nz*nhx))

  
  call off2ang2d_sd_init(z0, hx0, dz, dhx, nz, nhx, axmin, axmax, style)
  do imx=1, nmx
     data = 0.
     call sreed("in", data, 4*nz*nhx)
     stat = off2ang2d_sd_oper(T, F, modl, data)
  !call solver_tiny(m=modl, d=data, Fop=off2ang_sd_oper, stepper=cgstep, niter=50)
     call srite("out", modl, 4*nz*nhx)
  end do
  call to_history("n1", nz) ; call to_history("d1", dz) ; call to_history("o1", z0)
  call to_history("n2", nax); call to_history("d2", dax); call to_history("o2", axmin) 
  call to_history("n3", nmx); call to_history("d3", dmx); call to_history("o3", mx0)


end program