program multibunch

   use bmad
     implicit none

      type (lat_struct) lat
      type (ele_struct) ele 
      type (coord_struct), allocatable :: orbit(:)
       integer unit, number, m, i, ix, nargs,j
       integer ix_cache/0/
       integer lun
       integer nturns
       integer ios

      character*140 lat_file
      character*120 line
      logical error/.false./

      real(rp) offsetvec(1:6)

      namelist /run_def/nturns, offsetvec

  !this stuff is to ready an input file from the command line
      nargs = cesr_iargc()
      if(nargs == 1)then
         call cesr_getarg(1,lat_file)
         print *, 'Using ', trim(lat_file)
       else
        lat_file = 'bmad.'
        print '(a,$)',' Lattice file name ? (default= synch.lat) '
        read(5,'(a)') line
        call string_trim(line, line, ix)
        lat_file = line
        if(ix == 0) lat_file = 'synch.lat'
        print *, ' lat_file = ', lat_file
      endif

     lun = lunget()
     open(unit=lun, file='run_def.in', STATUS ='old')
     read(lun, nml=run_def, IOSTAT=ios)
     close(unit=lun)
     write(6,nml=run_def)

       call bmad_parser (lat_file, lat)     !parse lattice file
       call reallocate_coord (orbit, lat%n_ele_track)  !allocate memory for orbit
       call lat_make_mat6(lat, -1)               !compute transfer matrics for all of the elements in the lattice
       call twiss_at_start(lat)          ! compute twiss parameters at the starting point
       call twiss_propagate_all(lat)     ! propagate twiss parameters to all elements in the lattice


       call closed_orbit_calc(lat,orbit,6)     !compute the closed orbit

       orbit(0)%vec = orbit(0)%vec + offsetvec

     lun = lunget()
     open(unit=lun, file='turns.dat')

     write(lun,'(a10,6a12)')'turn','x','px','y','py','z','delta'
     do i=1,nturns       ! compute trajectory offset from the closed orbit for <nturn> turns
      call track_all(lat,orbit)
      orbit(0)%vec = orbit(lat%n_ele_track)%vec
      write(lun,'(i10,6es12.4)')i, orbit(lat%n_ele_track)%vec(1:6)
     end do

    close(unit=lun)

    call writeinfo(lat)    

end