subroutine first_turn(lat,spin_tracking_on, from_orbit, circumference, make_movie)
         use bmad
         use parameters_bmad
         use muon_interface, dummy => first_turn
         implicit none
         type (lat_struct) lat
         type (ele_struct) ele_at_s, ele_offset
         type (coord_struct), allocatable :: from_orbit(:), to_orbit(:)
         type (coord_struct) orb_at_s, co
         type (em_field_struct) field
         real(rp) circumference
         real(rp) Q(3)
         real(rp) s_end
         real(rp) s_start
         real(rp) s_rel
         real(rp) s, s_eff,c, step/0.01/
         real(rp) s_ele
         logical spin_tracking_on, make_movie
         logical local_ref/.false./,calcd/.false./
         logical err_flag
         logical eflag
         logical local_ref_frame/.false./
         integer j,i, ix_end
         integer ix_ele
         integer lun
         integer nbranch/1/
         integer ie_from, track_state
  
         i=1
                do j=1,lat%branch(0)%n_ele_track
                 if(spin_tracking_on) Q = from_orbit(j)%spin
                 call write_single_particle_by_element(i,lat%branch(0)%ele(j)%name,from_orbit(j), lat%branch(0)%ele(j)%s/circumference, &
                                                            ((i-1)*lat%branch(0)%ele(lat%branch(0)%n_ele_track)%s +lat%branch(0)%ele(j)%s)/circumference, Q)
                 call lat_make_mat6 (lat, j,from_orbit, ix_branch=0, err_flag=err_flag)
                end do

                call twiss_propagate_all(lat, ix_branch=0)
                do j=1,lat%branch(0)%n_ele_track
                  s_start=0.
                  if(make_movie)call step_floor_around_ring(lat,j,0,from_orbit(j-1),lat%branch(0)%ele(j-1)%floor%r, s_start)
                end do
 
                s_end=0.00001
                orb_at_s%vec=from_orbit(0)%vec
                ix_end=1
                lun=lunget()
                open(unit=lun, file=trim(directory)//'/'//'inj_traj_bfield.dat')
                write(lun,'(a,es12.4)')' inflector_field = ', inflector_field
                write(lun,'(11a12)')'s', 'x','xp','y','yp','z','dp','Bx','By','Bz','state'
                call write_injection_line_trajectory(s_end, orb_at_s%vec, from_orbit(1)%vec,lat%ele(1)%name, ix_end,lat%ele(1)%floor%r,orb_at_s%state)
                do while(s_end < lat%ele(lat%branch(0)%n_ele_track)%s +0.05) !the 0.05 is to make sure we overlap the very end
                   if(s_end > lat%ele(lat%branch(0)%n_ele_track)%s) s_end = lat%ele(lat%branch(0)%n_ele_track)%s !track to the very end

                   call twiss_and_track_at_s(lat, s_end, ele_at_s, from_orbit, orb_at_s, ix_branch=0,compute_floor_coords=.true.)          
                   print '(es12.4,1x,i10,1x,a,1x,es12.4,1x,a)', s_end, lat%branch(0)%n_ele_track,lat%ele(lat%branch(0)%n_ele_track)%name, lat%ele(lat%branch(0)%n_ele_track)%s, ele_at_s%name
                   ix_end = element_at_s (lat, s_end, .true.)
                   s_rel = s_end - lat%ele(ix_end-1)%s
                   call em_field_calc(lat%ele(ix_end), lat%param, s_rel, orb_at_s, local_ref, field,calcd, err_flag)
                   write(lun,'(es12.4,1x,6es12.4,3es12.4,1x,i12)')s_end, orb_at_s%vec, field%B(1:3), orb_at_s%state
                   call write_injection_line_trajectory(s_end, orb_at_s%vec, from_orbit(ix_end)%vec,lat%ele(ix_end)%name, ix_end, ele_at_s%floor%r,orb_at_s%state) !end_orb%vec
                   call write_injection_line_twiss(ele_at_s)

                   s_end = s_end + 0.01
                end do

! add this next bit 6/23/20 to propagate twiss through injection line and into ring
                nbranch = 1
                lat%branch(nbranch)%ele(0)%a%beta = ele_at_s%a%beta
                lat% branch(nbranch)%ele(0)%a%alpha = ele_at_s%a%alpha
                lat%branch(nbranch)%ele(0)%b%beta = ele_at_s%b%beta
                lat% branch(nbranch)%ele(0)%b%alpha = ele_at_s%b%alpha
                lat%branch(nbranch)%ele(0)%x%eta = ele_at_s%x%eta
                lat%branch(nbranch)%ele(0)%x%etap = ele_at_s%x%etap
                ie_from = lat%branch(nbranch)%ix_from_ele
                call twiss_propagate_all(lat,ix_branch = nbranch)
                !                from_orbit(0) = orb_at_s
                call reallocate_coord(to_orbit, lat%branch(nbranch)%n_ele_max)
                to_orbit(0) = from_orbit(ie_from)
                call track_all(lat, to_orbit, nbranch, track_state, err_flag)
                s=0.
                do while(s <= lat%branch(nbranch)%ele(lat%branch(nbranch)%n_ele_track)%s)
                   ix_ele= element_at_s(lat,s,choose_max=.true., ix_branch=nbranch, err_flag=eflag, s_eff=s_eff, position=co)
                   call twiss_and_track_at_s(lat,s,ele_at_s, to_orbit,  orb_at_s=orb_at_s, ix_branch=nbranch, use_last=.false.)
                   s_ele = lat%branch(nbranch)%ele(ix_ele)%s_start
                   call em_field_calc(lat%branch(nbranch)%ele(ix_ele), lat%param,s-s_ele,orb_at_s, local_ref_frame, field, err_flag=err_flag)
                   write(lun,'(es12.4,1x,6es12.4,3es12.4,1x,i12)')s_end, orb_at_s%vec, field%B(1:3), orb_at_s%state
                   ele_offset = ele_at_s
                   ele_offset%s =ele_at_s%s + s_end
                   call write_injection_line_twiss(ele_offset)           
!                   write(lun,'(i10,16es12.4,1x,a16,1x,3es12.4)')ix_ele,s+(i-1)*c,s_eff, ele_at_s%a%beta, ele_at_s%a%alpha, ele_at_s%b%beta, ele_at_s%b%alpha, ele_at_s%x%eta,&
!                                          ele_at_s%x%etap, ele_at_s%y%eta, ele_at_s%y%etap, lat%branch(nbranch)%ele(ix_ele)%a%beta,&
!                                          lat%branch(nbranch)%ele(ix_ele)%b%beta,orb_at_s%vec(1),orb_at_s%vec(3), &
!                                          closed_orb(ix_ele)%vec(1), closed_orb(ix_ele)%vec(3),lat%branch(nbranch)%ele(ix_ele)%name, field%b
                   s= s + step
               end do


                close(lun)

         return
         end