program read_vert_width_fit_params
  use sim_utils
  implicit none
  integer lun,i
  integer lun1
  integer reason
  integer j

  character*2 word
  character*4 word_e
  character*200 string
  character*20 try
  character(len=20) str
  real(rp) a(24),b(24),c(24),d(24),e(24)
  real(rp) f,x
  real(rp) location(24)
!  real(rp) location_s(24)  /4.8503,   8.0776,  1.0927E+01,  1.6022E+01,   1.9249E+01,   2.5083E+01,   2.7193E+01,   3.0421E+01,   3.3338E+01,  3.8365E+01,    4.1592E+01,  4.4686E+01/
  real(rp)location_s(24)/1.37,2.74,4.35,4.85,8.08,8.43,10.93,13.42,13.91,15.53,16.02,19.25,22.17,25.08,26.70,27.19,30.42,33.34,36.26,37.87,38.12,38.37,41.59,44.69/








  real(rp) g30(24), g150(24), g300(24)




  
  lun=lunget()
  open(unit=lun, file = 'plot_vertical_width_24.dat')
  do while(.true.)
     read(lun,'(a)', IOSTAT=reason)string
!     print '(a)', string
  if(reason <0) exit
  do i=1,24
     if(i<10)then
        write(word(1:1),'(i1)')i
        word_e = word(1:1)//' ='
     endif
     
     if(i>=10)then
        write(word(1:2),'(i2)')i
        word_e = word(1:2)//' ='
     endif
     try = str(i)
!     print '(a,a)',' try = ', try
!     print '(a10,a)',' word_e = ', word_e
     if (index(string, word_e) /=0 .and. string(1:1) == 'a')then
        read(string(6:),*)a(i)
     endif
     if (index(string, word_e) /=0 .and. string(1:1) == 'b') read(string(6:),*)b(i)
     if (index(string, word_e) /=0 .and. string(1:1) == 'c') read(string(6:),*)c(i)
     if (index(string, word_e) /=0 .and. string(1:1) == 'd') read(string(6:),*)d(i)
     if (index(string, word_e) /=0 .and. string(1:1) == 'e') read(string(6:),*)e(i)
  end do
end do
close(lun)
lun=lunget()
open(unit=lun, file = 'fitted_curves.dat')
lun1= lunget()
open(unit=lun1, file = 'vert_width_30_150us.dat')
write(lun1,'(3x,2a10, 5a12)') 'i','j','30','150','300','location','location_s'
print '(5a12)','a','b','c','d','e'
do i=1,24
   print '(5es12.4)',a(i),b(i),c(i),d(i),e(i)
   write(lun,'(/,a1,a12,1x,i2,/)')'#',  ' location = ', i
   do j=30,150
      x=float(j)
!      f = a(i) + b(i)*x + c(i)*x**2 + d(i)*x**3 + e(i)*x**4
     f = a(i)*exp(x/b(i)) + e(i)*exp(x/c(i)) +d(i)
      write(lun,'(i10,2es12.4)')i,x,f
      
   end do
   x=30.
   j=i
   g30(j)=a(j) + b(j)*x + c(j)*x**2 + d(j)*x**3 + e(j)*x**4
   g30(j)=a(j) * exp(x/b(j)) + e(j)*exp(x/c(j)) + d(j)
   
   x=150.
   g150(j)=a(j) + b(j)*x + c(j)*x**2 + d(j)*x**3 + e(j)*x**4
   g150(j)=a(j) * exp(x/b(j)) + e(j)*exp(x/c(j)) + d(j)   

   x=300.
   g300(j)=a(j) + b(j)*x + c(j)*x**2 + d(j)*x**3 + e(j)*x**4
   g300(j)=a(j) * exp(x/b(j)) + e(j)*exp(x/c(j)) + d(j)   
   !f1(x) = a1*exp(x/b1) + e1*exp(x/c1) +d1
   
   write(lun1, '(2i10,1x,5es12.4)')i,j,g30(j), g150(j), g300(j), location_s(j)/44.686*360., location_s(i) 
end do
print '(a)',' write "fitted_curves.dat"'
print '(a)', 'write "vert_width_30_150us.dat"'
end program read_vert_width_fit_params

character(len=20) function str(k)
!   "Convert an integer to string."
    integer, intent(in) :: k
    write (str, *) k
    str = adjustl(str)
end function str