*******************************************************************************

	program mpnm 

*******************************************************************************
c
c	This program takes DMol3 formatted  output and uses the data to
c       make output files that contain normal mode position shifts	
c
c	To compile non-graphics programs use : 
c             f77 $(OPT) mpnm.f -o mpnm 
c
c input

	character*72    atomform
        character*80    ag98form
        character*80    fg98form
        character*80    mg98form
        character*80    massform
	character*80    aform 
	character*80    record, rec
	character*80    recd(300)
	character*120   modeform
	character*120   record1
	character*80    freqform
        character*40	fname,gname
        character*40	hname,iname
        character*40	jname,kname,qname
        character*40	pname,rname
        character*10	oname
        character*4	root
        character*6	rname1
        character*7	rname2
        character*8	rname3
        character*1     tchar1,tempchar1
        character*2     tchar2,tempchar2
        character*3     tchar3,tempchar3
        character*2     AtomType(100) 
	real*8 		remain
	real*8 		rm(500)
	real*8 		xyz(3,500)
	real*8 		xd(500,500),yd(500,500)
	real*8 		xdn(500,500),ydn(500,500)
	real*8 		zd(500,500)
	real*8 		zdn(500,500)
	real*8          dd(3,500)
	real*8          freq(500)
	real*8          junk1(500)
	real*8          junk2(500)
	real*8          at(500)
        real*8          dev(500)
        real*8          devmax
        real*8          devinc
	integer  	ichoice
	integer  	isel(300)
	integer  	nsel,nn
	integer  	iatom,im,i,j,ierror
	integer  	ifreq,ii,il,ir
        integer         idmol	
	integer		iunit 
	integer		imass 
	integer		ipbc 
	integer		ilinear, iiso
	integer		ndev 
	integer	        iscale	
	integer		inum,icheck,itype
	integer  	iAtomType(100)
	integer         modeflag, nfreqflag	
	integer         iproj , ibeg
c
c 
        atomform = "(a2,11x,f16.12,4x,f16.12,4x,f16.12)"
        modeform = "(9x,9(f13.10))"
        freqform = "(3x,i2,4x,f8.6,5x,f7.1,6x,f6.2,39x)"
30      format(a2,11x,f16.12,4x,f16.12,4x,f16.12)
c
	 iunit=3
         print *, 'Enter 0 for G98 and 1 for Dmol3. '
         read(5,*) idmol 
         if(idmol.eq.1) then	
         print *, 'Enter outmol filename: '
         else if (idmol.eq.0) then
         print *, 'Enter G98 .out filename: '
         endif
         read(5,'(a40)') fname
c
c
         imass=0
         ipbc=0
         ilinear=0
         iiso=0
c
c print *,'Enter 1 for PBC=ON, Enter 0 for PBC=OFF: '
c        read(5,*) ipbc
c
         print *,'Enter 1 for linear, Enter 0 for non-linear: '
         read(5,*) ilinear
c
	 print *,'Enter number of atoms: '
         read(5,*) iatom
c
         if(ipbc.eq.0) then
           if(ilinear.eq.1) then
           ifreq=3*iatom-5
           else if(linear.eq.0) then 
           ifreq=3*iatom-6
           else
           print *,'Error> Linear flag is not set'
           stop
           endif
         endif
         if(ipbc.eq.1) then
          ifreq=3*iatom-3
        endif
c
         iscale=0
c
c iscale = 0 means mass-weighted coordinates
c
	 print *,'Enter maximum displacement :' 
         read(5,*) devmax
         devmax=abs(devmax)
	 print *,'Enter number of steps: ' 
         read(5,*) ndev
         devinc=2*devmax/real(ndev-1)
c         print *,devinc
        do i=1,ndev
          dev(i)=-devmax+(i-1)*devinc
        enddo
c
        modeflag=0
        nfreqflag=0
c         
c       Create a file that contains the optimized geometry from
c       the Gaussian output file immediately followed by the eigenvector
c       projections 
c
c        Read in Cartesian coordinates
c
        if(idmol.eq.1) then
        iproj=0
	 print *,'Enter 1 if vibration project is turned off. '
         read(5,*) iproj
        ibeg=0
        if(iproj.eq.1) then
          print *,'Vibration_project was turned off'
          print *,'Translations and rotations are included'
          print *,'in the projected eigenvectors.'
          if(ilinear.eq.0) then
            ifreq=ifreq+6
          else if(ilinear.eq.1) then
            ifreq=ifreq+5
          endif
        endif
        print *,'Enter 1 to include all eigenvectors.'
        read (5,*) ichoice
        if(ichoice.eq.0.and.ilinear.eq.0) ibeg=6
        if(ichoice.eq.0.and.ilinear.eq.1) ibeg=5
c
        print *,'ibeg = ',ibeg
        print *,'ifreq = ',ifreq
c         
         open(unit=3,name=fname,err=999)
        do 900 j=1,50000
	  read(iunit,'(a80)',end=999,err=999) record
        if(record.eq.'$coordinates') then
	 do i=1,iatom
	  read(iunit,'(a80)',end=999,err=999) record
          read (record,atomform,err=999) AtomType(i),xyz(1,i),
     1   xyz(2,i),xyz(3,i)
           if(AtomType(i).eq.'H') iAtomType(i)=1
           if(AtomType(i).eq.'B') iAtomType(i)=3
           if(AtomType(i).eq.'C') iAtomType(i)=6
           if(AtomType(i).eq.'N') iAtomType(i)=7
           if(AtomType(i).eq.'O') iAtomType(i)=8
           if(AtomType(i).eq.'F') iAtomType(i)=9
           if(AtomType(i).eq.'S') iAtomType(i)=16
           if(AtomType(i).eq.'Cl') iAtomType(i)=17
           if(AtomType(i).eq.'I') iAtomType(i)=53
           if(AtomType(i).eq.'Fe') iAtomType(i)=26
           if(AtomType(i).eq.'Ni') iAtomType(i)=28
           if(AtomType(i).eq.'Cu') iAtomType(i)=27
           if(AtomType(i).eq.'Ag') iAtomType(i)=47
           if(AtomType(i).eq.'Au') iAtomType(i)=79
           if(AtomType(i).eq.'H') at(i)=1.00794
           if(AtomType(i).eq.'B') at(i)=10.811
           if(AtomType(i).eq.'C') at(i)=12.0107
           if(AtomType(i).eq.'N') at(i)=14.00674
           if(AtomType(i).eq.'O') at(i)=15.9994
           if(AtomType(i).eq.'F') at(i)=18.9984
           if(AtomType(i).eq.'S') at(i)=32.066
           if(AtomType(i).eq.'Cl') at(i)=35.4527
           if(AtomType(i).eq.'I') at(i)=126.904
           if(AtomType(i).eq.'Fe') at(i)=55.845
           if(AtomType(i).eq.'Ni') at(i)=58.70
           if(AtomType(i).eq.'Cu') at(i)=63.546
           if(AtomType(i).eq.'Ag') at(i)=107.868
           if(AtomType(i).eq.'Au') at(i)=196.9665
           write(6,20) AtomType(i),xyz(1,i),xyz(2,i),xyz(3,i)
         enddo
        endif
c
        if(record.eq.'  mode     au_amu        cm-1      km/mol') then
          nfreqflag=1
        endif
        if(record.eq.'     vibrational frequencies') then
          read(iunit,'(a80)',end=999,err=999) record
          nfreqflag=2
        endif
        if(nfreqflag.gt.0.and.nfreqflag.le.2) then
	 do i=1,ifreq
	   read(iunit,'(a80)',end=999,err=999) record
        if(record.eq.'  Frequencies (cm-1) and normal modes') then
           goto 991
        endif
           write(6,'(a80)') record
          read(record,freqform,err=999) im,junk1(i),freq(i),junk2(i)
         enddo
         nfreqflag=3
        endif 
c
        if(record.eq.'  Frequencies (cm-1) and normal modes') then
         modeflag=1
        if(ifreq.ge.9) then
         il=ifreq/9
         remain=real(ifreq)/9.0-real(il)
         print *,'ir = ',ir
         if(remain.ne.0) il=il+1
       else
         il=1
         ir=0
       endif
         do 1100 ii=1,il
	  read(iunit,'(a80)',end=999,err=999) record
	  read(iunit,'(a80)',end=999,err=999) record
          i=(ii-1)*9+1
          do 1050 im=1,iatom
	  read(iunit,'(a120)',end=999,err=999) record1
	 read(record1,modeform,err=999) xd(i,im),xd(i+1,im),
     1 xd(i+2,im),xd(i+3,im),xd(i+4,im),xd(i+5,im),
     2 xd(i+6,im),xd(i+7,im),xd(i+8,im)
c       write(6,41) xd(i,im),xd(i+1,im),
c     1 xd(i+2,im),xd(i+3,im),xd(i+4,im),xd(i+5,im),
c     2 xd(i+6,im),xd(i+7,im),xd(i+8,im)
	  read(iunit,'(a120)',end=999,err=999) record1
	 read(record1,modeform,err=999) yd(i,im),yd(i+1,im),
     1 yd(i+2,im),yd(i+3,im),yd(i+4,im),yd(i+5,im),
     2 yd(i+6,im),yd(i+7,im),yd(i+8,im)
c       write(6,42) yd(i,im),yd(i+1,im),
c     1 yd(i+2,im),yd(i+3,im),yd(i+4,im),yd(i+5,im),
c     2 yd(i+6,im),yd(i+7,im),yd(i+8,im)
	  read(iunit,'(a120)',end=999,err=999) record1
	 read(record1,modeform,err=999) zd(i,im),zd(i+1,im),
     1 zd(i+2,im),zd(i+3,im),zd(i+4,im),zd(i+5,im),
     2 zd(i+6,im),zd(i+7,im),zd(i+8,im)
c       write(6,43) zd(i,im),zd(i+1,im),
c     1 zd(i+2,im),zd(i+3,im),zd(i+4,im),zd(i+5,im),
c     2 zd(i+6,im),zd(i+7,im),zd(i+8,im)
1050     continue
	  read(iunit,'(a80)',end=999,err=999) record
	  read(iunit,'(a80)',end=999,err=999) record
1100    continue
        goto 2000
        endif
900     continue
c
2000    continue 
        if(modeflag.eq.0) goto 992
        if(nfreqflag.eq.0) goto 993
        close(3)
c
c G98 calculation
c
        else if(idmol.eq.0) then
        ag98form = '(14x,i2,20x,3(f9.6,3x),8x)'
        fg98form = '(16x,3(f9.4,14x))'
        mg98form = '(11x,3(f7.3),2x,3(f7.3),2x,3(f7.3))'
        massform = '(38x,f9.5)'
        iproj=0
        icoord=0
        iopt=0
        do 915 j=1,50000
	  read(iunit,'(a80)',end=1999,err=999) record
          if(record(1:24).eq.' Optimization completed.') iopt=j
          if(record(26:33).eq.'Standard') icoord=j
915     continue
1999    continue
        print *,' iopt = ',iopt,' icoord = ',icoord
        close(3)
        open(unit=3,name=fname,err=999)
        if(iopt.gt.icoord) then
          istart=icoord
        else if(iopt.lt.icoord) then
          istart=iopt
        endif
          do 1915 j=1,istart
	  read(iunit,'(a80)',end=999,err=999) record
1915      continue
          print *, record
        if(istart.eq.iopt) then
          do 1917 j=1,1001
	  read(iunit,'(a80)',end=999,err=999) record
          if(record(26:33).eq.'Standard') goto 701
1917      continue
        endif
701     continue
c
        do 1919 j=1,4
	  read(iunit,'(a80)',end=999,err=999) record
1919    continue
         do i=1,iatom
	  read(iunit,'(a80)',end=999,err=999) record
          print *, record
          read (record,ag98form,err=999) iAtomType(i),xyz(1,i),
     1   xyz(2,i),xyz(3,i)
           if(iAtomtype(i).eq.1) AtomType(i)='H'
           if(iAtomtype(i).eq.3) AtomType(i)='B'
           if(iAtomtype(i).eq.6) AtomType(i)='C'
           if(iAtomtype(i).eq.7) AtomType(i)='N'
           if(iAtomtype(i).eq.8) AtomType(i)='O'
           if(iAtomtype(i).eq.26) AtomType(i)='Fe'
           if(iAtomType(i).eq.1) at(i)=1.00794
           if(iAtomType(i).eq.3) at(i)=10.811
           if(iAtomType(i).eq.6) at(i)=12.0107
           if(iAtomType(i).eq.7) at(i)=14.00674
           if(iAtomType(i).eq.8) at(i)=15.9994
           if(iAtomType(i).eq.26) at(i)=55.845
           write(6,20) AtomType(i),xyz(1,i),xyz(2,i),xyz(3,i)
         enddo
c
         close(3)
         open(unit=3,name=fname,err=999)
         do 905 j=1,50000
	  read(iunit,'(a80)',end=999,err=999) record
        if(record(1:13).eq.' Frequencies ') then
         nfreqflag=1
         do k=1,ifreq/3
          i=(k-1)*3+1     
          write(6,'(a80)') record
          read(record,fg98form,err=999) freq(i),freq(i+1),freq(i+2)
c
	  read(iunit,'(a80)',end=999,err=999) record
	  read(iunit,'(a80)',end=999,err=999) record
	  read(iunit,'(a80)',end=999,err=999) record
          read(iunit,'(a80)',end=999,err=999) record
	  read(iunit,'(a80)',end=999,err=999) record
	  read(iunit,'(a80)',end=999,err=999) record
          do 1055 im=1,iatom
          read(iunit,'(a80)',end=999,err=999) record
          read(record,mg98form,err=999) xd(i,im),yd(i,im),
     1 zd(i,im),xd(i+1,im),yd(i+1,im),zd(i+1,im),
     2 xd(i+2,im),yd(i+2,im),zd(i+2,im)
          write(6,41) xd(i,im),yd(i,im),
     1 zd(i,im),xd(i+1,im),yd(i+1,im),zd(i+1,im),
     2 xd(i+2,im),yd(i+2,im),zd(i+2,im)
1055     continue
	  read(iunit,'(a80)',end=999,err=999) record
	  read(iunit,'(a80)',end=999,err=999) record
	  read(iunit,'(a80)',end=999,err=999) record
        enddo
         endif
        if(nfreqflag.eq.1) goto 2993
905      continue
        if(nfreqflag.eq.0) goto 993
2993     continue
        close(3)
        if(iiso.eq.1) then
        open(unit=3,name=fname,err=999)
        do 907 j=1,50000
	  read(iunit,'(a80)',end=999,err=999) record
         if(record(1:6).eq.' Atom ') then
           imass=1
           do 919 jj=1,iatom 
           read(record,massform,err=999) at(jj)
	   read(iunit,'(a80)',end=999,err=999) record
           print *, at(jj)
919      continue
         endif
         if(imass.eq.1) goto 3993
907      continue
3993     continue
         endif
        endif
        if(idmol.eq.1) then
        print *,'Enter car filename: '
        else if(idmol.eq.0) then
        print *,'Enter G98 input filename: '
        endif
        read(5,'(a40)') hname
        open(unit=7,name=hname,err=999)

         aform = "(a5,3(2x,f13.9),a29)"
c
c DMol3 calculation
c
        if(idmol.eq.1) then
         read(7,'(a80)',end=997,err=995) recd(1)
         read(7,'(a80)',end=997,err=995) recd(2)
         read(7,'(a80)',end=997,err=995) recd(3)
         read(7,'(a80)',end=997,err=995) recd(4)
         if(ipbc.eq.1) then
         read(7,'(a80)',end=997,err=995) recd(5)
         endif

         do 200 i=1,iatom
           read(7,'(a80)',end=997,err=995) recd(i+4+ipbc)
           write(6,'(a80)') recd(i+4+ipbc)
           read (recd(i+4+ipbc),aform,err=996) left,xyz(1,i),
     $ xyz(2,i),xyz(3,i),right
200      continue
c
c G98 calculation
c
         else if (idmol.eq.0) then
c
c read header that consists of seven lines
c
         do 195 im=1,7 
         read(7,'(a80)',end=1997,err=995) recd(im)
         write(6,'(a80)') recd(im)
195      continue
         icnt=7
c
c read atomic coordinates from G98 input file
c
         icnt=icnt+iatom
         do 205 im=1,iatom
           read(7,'(a80)',end=1997,err=995) recd(im+icnt)
           write(6,'(a80)') recd(im+icnt)
c           read (recd(im+icnt),aform,err=996) left,xyz(1,i),
c     $ xyz(2,i),xyz(3,i),right
205      continue
         read(7,'(a80)',end=1997,err=995) recd(1+icnt)
         write(6,'(a80)') recd(1+icnt)
         icnt=icnt+1
         read(7,'(a80)',end=1997,err=995) recd(2+icnt)
         write(6,'(a80)') recd(2+icnt)
         icnt=icnt+1
         do 207 im=1,iatom
           read(7,'(a80)',end=1997,err=995) recd(im+icnt)
           write(6,'(a80)') recd(im+icnt)
           icnt=icnt+1
207      continue
1997     continue
         icnt=icnt-1
c         print *,'End of file detected in G98 input file'
         endif
         close(7)
c
         print *,'Enter rootname: '
         read(5,'(a4)') root 
c
	 print *,'Enter number of selected modes: '
         read(5,*) nsel
         if(nsel.ne.0) then 
         jname=root//'.modes'
         print *,'Reading file ',jname
         open(19,name=jname,form='formatted',err=999)
         do nn=1,nsel
           read(19,'(i4)') isel(nn)
         enddo
         close(19)
         endif 
c
c generate the mpnm file that contains general information 
c on the projections
c
        fname=root//'.mpnm'
        open(unit=17,file=fname,form='formatted',err=999)
          if(iscale.eq.0) write(17,330) 'mass-weighted              '
          if(iscale.eq.1) write(17,330) 'dimensionless mass-weighted'
          if(iscale.eq.2) write(17,330) 'normal coordinates only    '
          if(iproj.eq.1) then
            write(17,340) 'off.'
            if(ilinear.eq.1) then
              write(17,'(a49)') 'The first 5 files are translations and
     $ rotations.'
            else if(ilinear.eq.0) then
              write(17,'(a49)') 'The first 6 files are translations and
     $ rotations.'
             endif
          else if(iproj.eq.0) then 
            write(17,340) 'on. '
          endif 
            write(17,'(a13)') 'Displacements'
          do k=1,ndev
           write(17,67) dev(k)
          enddo
c
         rname=root//'.submit'
         open(unit=13,name=rname,err=999)
81       format('#!/bin/csh -f')
83       format('#@ error = ',a4,'.log')
85       format('#@ output= ',a4,'.log')
87       format('#@ notification = never')
89       format('#@ checkpoint = no')
90       format('#@ class = fewcpu')
91       format('#@ class = standard')
93       format('#@ job_cpu_limit = 400300')
95       format('#@ cpu_limit = 400000')
96       format('#@ node = 1')
97       format('#@ node = 2')
98       format('#@ total_tasks = 1')
99       format('#@ total_tasks = 4')
101      format('#@ job_type = parallel')
103      format('#@ network.MPI = css0,shared,US')
105      format('#@ queue')
107      format(a40)
         write(13,81)
         write(13,83) root
         write(13,85) root
         write(13,87)
         if(idmol.eq.1) write(13,91)
         if(idmol.eq.0) write(13,90)
         write(13,93)
         write(13,95)
         if(idmol.eq.1) then
         write(13,97)
         write(13,99)
         else if(idmol.eq.0) then
         write(13,96)
         write(13,98)
         endif
         write(13,101)
         write(13,103)
         write(13,105)
         open(unit=11,name=root,err=999)
         do k=1,ndev
           write(11,'(f8.5)') dev(k)
         enddo
         close(11)
c
c generate input files for multiple calculations
c
         if(idmol.eq.1) then
         oname=root//'.input'
         pname=root//'.inpgen'
         qname=root//'.prpgen'
         open(unit=15,name=pname,err=999)
         open(unit=9,name=qname,err=999)
77       format('cp ',a10,x,a20)
79       format('prepq ',a40)
78       format('g98q ',a40)
         do 2150 k=1,ndev
          call int2char(k,tchar1,tchar2,tchar3)
            if (k.le.9) then
               rname1 = root//tchar1//'_'
            else if(k.le.99) then
               rname2 = root//tchar2//'_'
            else if(k.le.999) then
               rname3 = root//tchar3//'_'
            endif
        do 1080 i=1,ifreq
          write(17,56) i
          do 1070 im=1,iatom
           xdn(i,im) = xd(i,im)/sqrt(at(im))*dev(k)/0.529177
           ydn(i,im) = yd(i,im)/sqrt(at(im))*dev(k)/0.529177
           zdn(i,im) = zd(i,im)/sqrt(at(im))*dev(k)/0.529177
           write(17,57) xdn(i,im),ydn(i,im),zdn(i,im)
1070    continue
1080    continue
	do i=1+ibeg,ifreq
          call int2char(i,tchar1,tchar2,tchar3)
            if (i.le.9) then
              if(k.le.9) then
               gname = rname1//tchar1//'.car'
               iname = rname1//tchar1//'.input'
               jname = rname1//tchar1
               kname = rname1//tchar1//'.job'
              else if(k.le.99) then
               gname = rname2//tchar1//'.car'
               iname = rname2//tchar1//'.input'
               jname = rname2//tchar1
               kname = rname2//tchar1//'.job'
              else if(k.le.999) then
               gname = rname3//tchar1//'.car'
               iname = rname3//tchar1//'.input'
               jname = rname3//tchar1
               kname = rname3//tchar1//'.job'
              endif
            else if(i.le.99) then
              if(k.le.9) then
               gname = rname1//tchar2//'.car'
               iname = rname1//tchar2//'.input'
               jname = rname1//tchar2
               kname = rname1//tchar2//'.job'
              else if(k.le.99) then
               gname = rname2//tchar2//'.car'
               iname = rname2//tchar2//'.input'
               jname = rname2//tchar2
               kname = rname2//tchar2//'.job'
              else if(k.le.999) then
               gname = rname3//tchar2//'.car'
               iname = rname3//tchar2//'.input'
               jname = rname3//tchar2
               kname = rname3//tchar2//'.job'
              endif
            else if(i.le.999) then
              if(k.le.9) then
               gname = rname1//tchar3//'.car'
               iname = rname1//tchar3//'.input'
               jname = rname1//tchar3
               kname = rname1//tchar3//'.job'
              else if(k.le.99) then
               gname = rname2//tchar3//'.car'
               iname = rname2//tchar3//'.input'
               jname = rname2//tchar3
               kname = rname2//tchar3//'.job'
              else if(k.le.999) then
               gname = rname3//tchar3//'.car'
               iname = rname3//tchar3//'.input'
               jname = rname3//tchar3
               kname = rname3//tchar3//'.job'
             endif
            endif
	  do j=1,iatom
	    dd(1,j)=xyz(1,j)-xdn(i,j)
	    dd(2,j)=xyz(2,j)-ydn(i,j)
	    dd(3,j)=xyz(3,j)-zdn(i,j)
c            print *, xd(i,j),yd(i,j),zd(i,j)
	  enddo
           write(15,77) oname,iname
           write(9,79) jname
           if(nsel.ne.0) then
             do nn=1,nsel
               if(i.eq.isel(nn)) then
                 write(13,107) kname
               endif
           enddo
           else if(nsel.eq.0) then
               write(13,107) kname
           endif
           junit=7
           open(unit=7,file=gname,form='formatted',err=999)
           do 50 j = 1,4
             write(junit,'(a80)') recd(j)
50         continue 
           if(ipbc.eq.1) then
             write(junit,'(a80)') recd(5)
           endif
           do 100 j = 1,iatom
             rec=recd(j+4+ipbc)
           write (junit,20) 
     $ rec(1:8),dd(1,j),dd(2,j),dd(3,j),rec(51:80)
100     continue
          do 150 j = 1,2
	    write(junit,'(a3)') 'end'
150       continue 
        close(7)
        enddo
2150    continue
        close(15)
        close(9)
        else if(idmol.eq.0) then
         oname=root//'.inp'
         qname=root//'.prpgen'
         open(unit=9,name=qname,err=999)
         do 2155 k=1,ndev
          call int2char(k,tchar1,tchar2,tchar3)
            if (k.le.9) then
               rname1 = root//tchar1//'_'
            else if(k.le.99) then
               rname2 = root//tchar2//'_'
            else if(k.le.999) then
               rname3 = root//tchar3//'_'
            endif
        do 1085 i=1,ibeg+ifreq
          write(17,56) i
          do 1075 im=1,iatom
           xdn(i,im) = xd(i,im)/sqrt(at(im))*dev(k)
           ydn(i,im) = yd(i,im)/sqrt(at(im))*dev(k)
           zdn(i,im) = zd(i,im)/sqrt(at(im))*dev(k)
           write(17,57) xdn(i,im),ydn(i,im),zdn(i,im)
1075    continue
1085    continue
	do i=1+ibeg,ifreq
          call int2char(i,tchar1,tchar2,tchar3)
            if (i.le.9) then
              if(k.le.9) then
               gname = rname1//tchar1//'.inp'
               jname = rname1//tchar1
               kname = rname1//tchar1//'.job'
              else if(k.le.99) then
               gname = rname2//tchar1//'.inp'
               jname = rname2//tchar1
               kname = rname2//tchar1//'.job'
              else if(k.le.999) then
               gname = rname3//tchar1//'.inp'
               jname = rname3//tchar1
               kname = rname3//tchar1//'.job'
              endif
            else if(i.le.99) then
              if(k.le.9) then
               gname = rname1//tchar2//'.inp'
               jname = rname1//tchar2
               kname = rname1//tchar2//'.job'
              else if(k.le.99) then
               gname = rname2//tchar2//'.inp'
               jname = rname2//tchar2
               kname = rname2//tchar2//'.job'
              else if(k.le.999) then
               gname = rname3//tchar2//'.inp'
               jname = rname3//tchar2
               kname = rname3//tchar2//'.job'
              endif
            else if(i.le.999) then
              if(k.le.9) then
               gname = rname1//tchar3//'.inp'
               jname = rname1//tchar3
               kname = rname1//tchar3//'.job'
              else if(k.le.99) then
               gname = rname2//tchar3//'.inp'
               jname = rname2//tchar3
               kname = rname2//tchar3//'.job'
              else if(k.le.999) then
               gname = rname3//tchar3//'.inp'
               jname = rname3//tchar3
               kname = rname3//tchar3//'.job'
             endif
            endif
	  do j=1,iatom
	    dd(1,j)=xyz(1,j)-xdn(i,j)
	    dd(2,j)=xyz(2,j)-ydn(i,j)
	    dd(3,j)=xyz(3,j)-zdn(i,j)
c            print *, dd(1,j),dd(2,j),dd(3,j)
	  enddo
           write(unit=9,78) jname
           if(nsel.ne.0) then
             do nn=1,nsel
               if(i.eq.isel(nn)) then
                 write(13,107) kname
               endif
           enddo
           else if(nsel.eq.0) then
               write(13,107) kname
           endif
           junit=7
           open(unit=7,file=gname,form='formatted',err=999)
           do 55 j = 1,7
             write(junit,'(a80)') recd(j)
55         continue 
           do 117 j = 1,iatom
           write (junit,28) AtomType(j),dd(1,j),dd(2,j),dd(3,j)
117     continue
           do 157 j = iatom+8,icnt
             write(junit,'(a80)') recd(j)
157         continue 
         write(junit,'(a3)') '   '
        close(7)
        enddo
2155    continue
        close(9)
        endif
        close(13)
        close(17)
c
330      format('The scaling option is ',a26)
340      format('Vibration projection is ',a4)
20       format(a7,f13.9,2x,f13.9,2x,f13.9,a30)
28      format(a2,3x,3(f9.6,3x))
41       format('x',9(x,f5.3))
42       format('y',9(x,f5.3))
43       format('z',9(x,f5.3))
56       format('Mode ',i2)
57       format(3(f12.3,x))
67       format(f8.3)
        goto 3001
995      continue
         print *,'Error in format'
         print *, record
         goto 999
991      continue
         print *,'Check number of atoms or PBC flag'
         print *,'Program expects too many frequencies'
         print *,'ifreq = ',ifreq
         goto 999
992      continue
         print *,'No normal modes found'
         goto 999
993      continue
         print *,'No frequencies found'
         goto 999
996      continue
         print *,'Error parsing format'
         print *, record
         print *,left, xyz(1,i), xyz(2,i), xyz(3,i),right 
         goto 999
997      continue
         print *,'End of file detected'
999      continue
        print *,'Error in file'
        goto 3002
3001    continue
        print *,'PROGRAM COMPLETE'
        print *,'DMOL INPUT FILES HAVE BEEN GENERATED'
3002    continue
	end

***************************************************************************

      Subroutine int2char(i,tempchar1,tempchar2,tempchar3)

***************************************************************************
!          To convert the integer data to character data

      Implicit None

      Integer i
      Character*1 tempchar1
      Character*2 tempchar2
      Character*3 tempchar3
      Character*4 record
 
      if(i.lt.10) then
        Write (record,'(i1)') i
        Read (record,'(a1)') tempchar1
      else if(i.lt.100) then
        Write (record,'(i2)') i
        Read (record,'(a2)') tempchar2
      else if(i.lt.1000) then
        Write (record,'(i3)') i
        Read (record,'(a3)') tempchar3
      endif

      Return
      End