!+ 
! Program wall_generator
!
! Parses lattice file and outputs wall.out file with points that
! can be linearly interpolated to make wall geometry
! Note: has header for particle 1, as if it were particle track data
!
! Modules Needed:
!   use bmad
!
! Input (command line)
!   [lattice file]
!
! Output
!   wall.dat
!-

program wall_generator

use bmad
use capillary_mod !For wall perpendicular


implicit none

type (lat_struct) :: lat
type (ele_struct), pointer :: ele
type (coord_struct) :: point
type (photon_coord_struct) :: p_orb
character(100) :: lat_name
integer :: ix_ele, section_id
integer, parameter :: wall_file = 1
real(rp) :: perp(3), dummy_radius

!------------------------------------------
!Get data file name
call getarg(1, lat_name)
print *, "Creating wall for lattice file: ", lat_name

!Parse lattice
call bmad_parser (lat_name, lat)

!-----------------------
!Open output file
open(wall_file, file = "wall.out")

!Write to file: header
write (wall_file, '(a, i8)'),   'Wall'
write (wall_file, '(8a18)'), 't', 'x', 'normal_x', 'y', 'normal_y', &
     's', 'normal_s'
write (wall_file, '(8a18)'), 's', 'm', '1', 'm', '1', &
     'm', '1'


!Make wall border points and write to file
!First set invariant coordinates
point%t = 0.0
point%vec(2) = 0.0
point%vec(4) = 0.0
point%vec(6) = 0.0
point%phase_x = 0.0
point%vec(3) = 0.0


!Set rest of coordinates
do ix_ele = 0, lat%n_ele_track

   ele => lat%ele(ix_ele)

   if (associated(ele%wall3d%section)) then

      !Move along wall cross sections and get coordinates at y = 0
      do section_id = 1, size(ele%wall3d%section)
         point%vec(1) = ele%wall3d%section(section_id)%v(1)%radius_x
         point%vec(5) = ele%wall3d%section(section_id)%s
         
		 !Get wall perpendicular
		 p_orb%orb%vec = point%vec
		 dummy_radius =  capillary_photon_d_radius (p_orb, ele, perp)
         point%vec(2) = perp(1)
         point%vec(4) = perp(2)
         !Shift s position to global coordinates
         point%vec(5) =    point%vec(5) + ele%s - ele%value(l$)
         point%vec(6) = perp(3)
                   
         !Write to file: wall border point
         write (wall_file, '(8es18.10)') point%t, point%vec, point%phase_x
      end do

   else

      !Get coordinates at beginning aperture
      point%vec(1) = ele%value(x1_limit$)
      point%vec(5) = ele%s - ele%value(l$)

      !Write to file: beginning point at aperture
      write (wall_file, '(8es18.10)') point%t, point%vec, point%phase_x

      !Get coordinates at end aperture
      point%vec(1) = ele%value(x1_limit$)
      point%vec(5) = ele%s

      !Write to file: end point at aperture
      write (wall_file, '(8es18.10)') point%t, point%vec, point%phase_x

   end if

end do


!Close data file
close(wall_file)
print *, "Written: wall.out"

end program