program test

use bmad
implicit none
integer i,j
real(rp) qx/0.121/, qz/0.116/
real(rp) xc(0:100000), xr(0:100000), xm(0:100000)
real(rp) y(0:100000)
real(rp) a,b,alpha, c, qzp
real(rp) gain/20./
real(rp) delta_phase, dphi

  alpha = 0.01
 print '(a,$)',' Gain ?'
 accept *,gain
 print '(a,$)', ' alpha ? '
 accept *,alpha
  a = 1.
  b=1.


  y(0) = 0.
  qzp=qz
  do i=1,100000
   xc(i) = a * sin(twopi*qx*i)
   xr(i) = b * cos(twopi*qzp*i + gain*y(i-1))

   delta_phase = twopi*i*(qx-qzp)- gain*y(i-1)
   dphi = gain*y(i-1)
   if(abs(dphi) > twopi)then
      dphi = dphi - sign(twopi, dphi)
      y(i-1) = dphi/gain
   endif
!   xr(i) = b * cos(twopi*(qzp + gain*y(i-1))*i)
   xm(i) = xc(i)*xr(i)

   call low_pass(alpha,xm(i), y(i), y(i-1))

!   qzp = qz - gain*y(i)/twopi/5
! 

   write(11, '(i,7e12.4)')i,xc(i), xr(i), xm(i), y(i), qzp,delta_phase, dphi
 end do

end