c
c  Driver program used in seisan for fpfit
c
c  J havskov, august 2010

      implicit none
      character*80 text    ! input text
      character*80 t_out   ! output text
      integer i

c
c   make parameter file
c
      call make_fpfit_par

c
c   make run file
c
      open(33,file='fpfit.run',status='unknown')
      write(33,'(a)') 'fps'
      write(33,'(a)') 'sto'
      write(33,'(a)') '   '
      close (33)
c
c   convert seisan print.out to fpfit h71 print file used for input to fpfit
c
      open(33,file='print.out',status='old')
      open(34,file='fpfit.dat',status='unknown')
      i=0
c
c   find header line
c
 1    continue
      read(33,'(a)',end=99) text
      if(text(4:7).eq.'date') goto 2
      goto 1
 2    continue
      read(33,'(a)') text
      write(34,'(a,a)')
     *'  DATE    ORIGIN   LATITUDE',
     *' LONGITUDE  DEPTH    MAG NO           RMS'
c
      t_out=' '
      t_out(1:21)=text(1:21)
      if(text(27:27).eq.'N') t_out(22:22)='n'
      if(text(27:27).eq.'S') t_out(22:22)='s'
      t_out(23:27)=text(22:26)
      t_out(29:31)=text(29:31)
      if(text(37:37).eq.'E') t_out(32:32)='e'
      if(text(37:37).eq.'W') t_out(32:32)='w'
      t_out(33:37)=text(32:36)
      t_out(39:44)=text(38:43)                    ! depth
      t_out(49:51)='1.0'                          ! fix to 1.0
      t_out(52:54)=text(46:48)
      t_out(64:68)=text(53:57)                    ! rms
      write(34,'(a)') t_out
c
c   blank line
c
      write(34,'(a)') ' '
c
c  text
c
      write(34,'(a,a)')'  STN  DIST  AZ',
     *' TOA PRMK HRMN  PSEC TPOBS              PRES  PWT'
c
c  search input phase lines
c
 3    continue
      read(33,'(a)',end=99) text
      if(text(2:4).eq.'stn') goto 4
      goto 3
c
 4    continue
      read(33,'(a)',end=99) text
      if(text.ne.' ') then
         t_out=' '
         t_out(1:15)=text(1:15)
         t_out(17:20)=text(18:20)               ! ain
         if(text(26:26).eq.'P') then            ! only use P-phases
            t_out(22:22)=text(26:26)
         else
            go to 4
         endif
         t_out(23:23)=text(32:32)               ! pol, check if PP etc ?
         t_out(26:35)=text(42:51)
         t_out(54:59)=text(66:71)
         t_out(62:65)=text(73:76)
         t_out(36:41)=text(52:57)
         write(34,'(a)') t_out
         i=i+1
         goto 4
      else
         goto 3        ! possibely add next event to composite solution
      endif 
c
c   end of file
c
 99   continue
      write(6,*)' Number of phases',i
      close(33)
      close(34)
c
c   run fpfit
c
      call systemc('fpfit<fpfit.run',15)

      stop
      end


      subroutine make_fpfit_par
c
c   make fpfit parmeter file, defualt name is fpfit.inp
c
      open(33,file='fpfit.inp',status='unknown')

      write(33,'(a)') "ttl   1 'none'"
      write(33,'(a)') "hyp 'fpfit.dat'"
      write(33,'(a)') "out 'fpfit.out'"
      write(33,'(a)') "sum 'fpfit.fps'"
      write(33,'(a)') "pol 'fpfit.pol'"
      write(33,'(a)') 'for   1'
      write(33,'(a)') 'mag   0.000'    
      write(33,'(a)') 'obs   6'
      write(33,'(a)') 'dis  0.1000E+06'
      write(33,'(a)') 'res   100.0'    
      write(33,'(a)') 'ain   0.000      180.0'    
      write(33,'(a)') 'amp   0'
      write(33,'(a)') 'bst   0'
      write(33,'(a)') 'fin   1'
      write(33,'(a)') 'rep   1'
      write(33,'(a)') 'cmp   0'
      write(33,'(a)') 'hdr  0.1000E-01 0.2000E-01 0.5000E-01 0.1000'    
      write(33,'(a)') 'mcr   1.000      1.000      1.000      1.000'   
      write(33,'(a)') 'dir   90.00      250.0      10.00      1.000'   
      write(33,'(a)') 'dip   0.000      90.00      10.00      1.000'    
      write(33,'(a)') 'rak  -180.0      180.0      20.00      1.000'    
      write(33,'(a)') "chn   1 'VHZ'"
      write(33,'(a)') "hds   1 '2'"
      close(33)
      return
      end