PROGRAM MAKEBONDSLONG
      IMPLICIT none
      integer*4 nt0
      parameter (nt0=23000)
      character*80 fn,name
      real*4 bonding(88),r(nt0,3),rb,x1,y1,z1,rt,x2,y2,z2,rfac
      parameter (rfac=1.15)
      integer*4 i,nat,j,nbt(nt0,7),iz(nt0),nb(nt0),iargc,
     1ti,tj
c       Du    H   He   Li   Be    B    C    N    O    F
c       Ne   Na   Mg   Al   Si    P    S   Cl   Ar    K   Ca   Sc   Ti
c        V   Cr   Mn   Fe   Co   Ni   Cu   Zn   Ga   Ge   As 
      data bonding/
     10.50,0.32,0.98,1.28,0.95,0.87,0.82,0.80,0.78,0.77,
     10.76,1.59,1.41,1.23,1.16,1.11,1.20,1.04,1.03,2.08,1.79,1.49,1.37,
     11.27,1.23,1.22,1.22,1.21,1.20,1.22,1.30,1.31,1.27,1.25,1.21,1.19,
     11.17,2.21,1.96,1.67,1.50,1.39,1.35,1.32,1.30,1.30,1.33,1.39,1.53,
     11.49,1.46,1.45,1.41,1.38,1.36,2.40,2.03,1.74,1.70,1.70,1.69,1.68,
     11.67,1.90,1.66,1.64,1.64,1.63,1.62,1.61,1.75,1.61,1.49,1.39,1.35,
     11.33,1.31,1.32,1.35,1.39,1.54,1.53,1.52,1.51,1.51,1.50,1.49,0.10/
      
c
      WRITE(6,3000)
3000  FORMAT(' Makes bond table to an MCM coordinate file',/,
     1       '   (longer bonds for MD snapshots)',/,/,
     1       '   (old file will be rewritten!!)',/,/,
     1       'Filename:') 
      if(iargc().gt.0)then
       call getarg(1,fn)
      else
       read(5,'(a)')fn
      endif
      open(2,file=fn)
      read(2,2000)name
2000  format(a80)
      read(2,*)nat
      if(nat.gt.nt0)call report('too many atoms')
      do 1 i=1,nat
1     read(2,*)iz(i),(r(i,j),j=1,3)
      write(6,*)nat,' atoms'

      do 4 i=1,nat
      nb(i)=0
      do 4 j=1,7
4     nbt(i,j)=0
c
      do 3 i=1,nat
      ti=iz(i)+1
      if(ti.lt.1.or.ti.gt.89)call report('type out of range')
      rb=bonding(ti)*rfac
      x1=r(i,1)
      y1=r(i,2)
      z1=r(i,3)
      do 3 j=i+1,nat
      tj=iz(j)+1
      if(ti.ne.2.or.tj.ne.2)then
       if(tj.lt.1.or.tj.gt.88)call report('type out of range')
       rt=(rb+bonding(tj)*rfac)**2
       if(i.eq.86.and.j.eq.87)then
        write(6,*)'bo i rfac rb',bonding(ti),rfac,rb
        write(6,*)'bo j rfac rt',bonding(tj),rfac,rt
        write(6,*)'d2',(x1-r(j,1))**2+(y1-r(j,2))**2+(z1-r(j,3))**2
       endif
       x2=(x1-r(j,1))**2
       if(x2.lt.rt)then 
        y2=(y1-r(j,2))**2+x2
        if(y2.lt.rt)then
         z2=(z1-r(j,3))**2+y2
         if(z2.lt.rt)then
          nb(i)=nb(i)+1
          nb(j)=nb(j)+1
          nbt(i,nb(i))=j
          nbt(j,nb(j))=i
         endif
        endif
       endif
      endif
3     continue
c
      rewind 2
      write(2,2000)name
      write(2,*)nat
      do 2 i=1,nat
2     write(2,20)iz(i),(r(i,j),j=1,3),(nbt(i,j),j=1,7),0.0d0
20    format(i2,3f15.8,7i7,f6.2)
      close(2)
      write(6,*)' File ',fn,' was rewritten.'
      stop
      end

      subroutine report(s)
      character*(*) s
      write(6,*)s
      stop
      end