program tscan
      implicit none
      integer*4 NW,i,NEQUI,NENS,NEPRO,i1,i2,i3,i4,
     1istart,iend,j,ist,ien
      real*8 xmin,xmax,TEMP,NFS,NPS,PRES,x,dx,b
      character*10 so,sop
      character*80 s80,dynamic
      logical lex

c
      dynamic='minimize'
      XMIN=0.0d0
      XMAX=0.0d0
      NW=50
    
      write(6,600)
600   format(' Makes the driving script for Tinker torsion scan',/,
     1       '       Input: HAM.OPT',/,
     1       '              KEY.TXT (the .key file torso)',/,
     1       '              ham.xyz (for MD)',/,
     1       '              ham.prm (for MD)',/,
     1       '      Output: goham',/)
      open(9,file='HAM.OPT')
1     read(9,900,end=99,err=99)so
900   format(a10)
      write(6,900)so
      if(so(1:4).eq.'DYNA')read(9,'(a)')dynamic
      if(so(1:4).eq.'XMIN')read(9,*)xmin
      if(so(1:4).eq.'XMAX')read(9,*)xmax
      if(so(1:2).eq.'NW'  )read(9,*)NW
      if(so(1:4).eq.'END')goto 99
      goto 1
99    close(9)
      write(6,*)'Option read'
      inquire(file='ham.xyz',exist=lex)
      if(.not.lex)call report('ham.xyz not found!')
      inquire(file='ham.prm',exist=lex)
      if(.not.lex)call report('ham.prm not found!')
      inquire(file='KEY.TXT',exist=lex)
      if(.not.lex)call report('KEY.TXT not found!')
c
c     List of options
c     ---------------
c     dynamic the dynamic command, with system path (e.g. /bin/dynamic)
c     xmin lower coordinate limit
c     xmax upper coordinate limit
c     NEQUI number of equilibrium steps
c     NEPRO number of production steps
c     NENS  ensemble type (2-NVT, 4-NPT)
c     NFS   MD integration step in fs
c     NPS   MD time between dumps in ps --- only one dump per simulation
c     NINT histogram integration step -- not used here
c     NW   number of weighted histograms
c     PRES pressure (for NPT)
c     TEMP temperature
c
c     KEY.TXT typically contains periodic boundary condition
c     definition (A-AXIS ...) and the restricted coordinate
c
      open(9,file='goscan')
      dx=(xmax-xmin)/dble(NW-1)
      x=xmin-dx
      write(9,*)'rm ham.xyz_* *.dyn ham.0* ham.1* ham.2* ham.3*
     1ham.4* ham.5* ham.6* ham.7* ham.8* ham.9*'
      do 2 i=1,NW
c
      write(so,'(i10)')i
      do 21 istart=1,10
21    if(so(istart:istart).ne.' ')goto 22
22    do 23 iend=10,1,-1
23    if(so(iend:iend).ne.' ')goto 24
24    x=x+dx
c
c     prepare key file
      open(91,file='ham.key.'//so(istart:iend))
3     open(92,file='KEY.TXT')
      read(92,980,end=33,err=33)s80
980   format(a80)
      if(s80(1:17).eq.'RESTRAIN-DISTANCE')then
       read(s80(18:80),*)i1,i2,b
       write(91,9191)i1,i2,b,x,x
9191   format('RESTRAIN-DISTANCE',2i5,3f12.4)
      else
       if(s80(1:16).eq.'RESTRAIN-TORSION')then
        read(s80(18:80),*)i1,i2,i3,i4,b
        write(91,9192)i1,i2,i3,i4,b,x,x
9192    format('RESTRAIN-TORSION',4i5,3f12.4)
       else
        write(91,980)s80
       endif
      endif
      goto 3
33    close(91)
      close(92)

      write(9,'(a)')'cp ham.key.'//so(istart:iend)//' ham.key'
c
      do 211 ist=1,80
211   if(dynamic(ist:ist).ne.' ')goto 221
221   do 233 ien=80,1,-1
233   if(dynamic(ien:ien).ne.' ')goto 241
241   write(9,9031)(dynamic(j:j),j=ist,ien)
9031  format(1a1,$)
      write(9,903)
903   format(' ham ',' 0.01')
2     write(9,*)'cp ham.xyz '//so(istart:iend)//'.xyz'

      close(9)
      call report('goscan written')
      end

      subroutine report(s)
      character*(*) s
      write(6,60)
      write(6,'(a)')s
      write(6,60)
60    format(70(1h*))
      stop
      end