Subroutine Track_with_rf_modulation

		
subroutine  Track_with_rf_modulation
use madx_ptc_module
use pointer_lattice
use gauss_dis
implicit none
type(layout), pointer :: ALS
type(fibre), pointer ::p
integer i,mf,j,n
type(normal_spin) normal
type(vecresonance) vecres
type(probe) Ray
type(probe_8) Ray_8
type(damapspin) m
integer pos
type(internal_state) state
real(dp) DC_ac,A_ac,theta_ac,x(6),deps

DC_ac=1.0d0
A_ac=1.d-3
theta_ac=0.01d0

als=>M_U%start

als%last=>als%start
als%lastpos=als%last%pos

call move_to(als,p,"BEND",pos)
write(6,*) pos,p%mag%name,p%mag%vorname

                if(.not.associated(P%MAG%DC_ac)) then
                   allocate(P%MAG%DC_ac)
                   allocate(P%MAG%A_ac)
                   allocate(P%MAG%theta_ac)

                   allocate(P%MAGP%DC_ac)
                   allocate(P%MAGP%A_ac)
                   allocate(P%MAGP%theta_ac)
                   CALL alloc(P%MAGP%DC_ac)
                   CALL alloc(P%MAGP%A_ac)
                   CALL alloc(P%MAGP%theta_ac)

                   P%MAG%DC_ac=DC_ac
                   P%MAG%A_ac=A_ac
                   P%MAG%theta_ac=theta_ac*twopi
                   P%MAGP%DC_ac=DC_ac
                   P%MAGP%A_ac=A_ac
                   P%MAGP%theta_ac=theta_ac*twopi
                   P%MAG%slow_ac=.true.
                   P%MAGP%slow_ac=.true.
                ENDIF
    state=default0+only_4d0
 CALL FIND_ORBIT_x(als,X,state,1.0e-5_dp,fibre1=1)
  Ray=x
  ray%ac%om=0.01d0
   ray%ac%x(1)=1.d0    ! fake phasor 
   ray%ac%x(2)=0.d0    !   fake phasor 
   
    state=state+modulation0
    write(6,*) "Tracking data in modulation.dat"
    call kanalnummer(mf,"modulation.dat")
   write(6,*) ray%x
   write(6,*) ray%ac%x
          WRITE(mf,'(4(1x,E15.8))') ray%x(1:4)
   do i=1,10000
    call TRACK_PROBE(als,ray,state, FIBRE1=1)
       WRITE(mf,'(4(1x,E15.8))') ray%x(1:4)
  enddo
  close(mf)
    write(6,*) "Normal form data in modulation.dat"
    call kanalnummer(mf,'modulation_result.txt')
   write(6,*) ray%x
   write(6,*) ray%ac%x
   
    CALL INIT(state,3,0)

   call alloc(m)
   call alloc(ray_8)
   call alloc(normal)
   call alloc(vecres)
     Ray=x
   ray%ac%om=0.001d0
   ray%ac%x(1)=1.d0
   ray%ac%x(2)=0.d0

   m=1
   ray_8=ray+m
    Write(mf,*) " Initial  Probe_8 Ray_8 "
      call print(ray_8,mf)
      call TRACK_PROBE(als,ray_8,state, FIBRE1=1)
   Write(mf,*) " Final Probe_8 Ray_8 "
   call print(ray_8,mf)
   m=ray_8
     Write(mf,*) " Final damapspin m"
    call print(m,mf)
    normal=m
   
      Write(mf,*) " Tunes  including  RF modulation"
      write(mf,*) normal%n%tune
      vecres=normal%n%normal%nonlinear
      Write(mf,*) "   "
      Write(mf,*) " Sine Part of Vector Field of Normal Form in resonance basis (Tunes) "
      Write(mf,*) " including  RF modulation "
      deps=1.d-8 
      call print(vecres%sin,mf,deps)


   close(mf)
   
   call kill(m)
   call kill(ray_8)
   call kill(normal)
   call kill(vecres)
    
end subroutine  Track_with_rf_modulation