PROGRAM PAR
C     extracts vibrational data (geometry, FF) from CPMD output
      IMPLICIT none
      integer*4 NT0,MX3,NAT,IX,IY,IA,K,ii,N
      PARAMETER (NT0=500,MX3=3*NT0)
      real*8 F(MX3,MX3),X(3,NT0)
      LOGICAL lex
      character*2 atsy(89),AT(NT0)
      data atsy/' H','He','Li','Be',' B',' C',' N',' O',' F','Ne',
     3'Na','Mg','Al','Si',' P',' S','Cl','Ar',
     4' K','Ca','Sc','Ti',' V','Cr','Mn','Fe','Co','Ni','Cu','Zn',
     4          'Ga','Ge','As','Se','Br','Kr',
     5'Rb','Sr',' Y','Zr','Nb','Mo','Tc','Ru','Rh','Pd','Ag','Cd',
     5          'In','Sn','Sb','Te',' I','Xe',
     6'Cs','Ba','La',
     6               'Ce','Pr','Nd','Pm','Sm','Eu','Gd','Tb','Dy','Ho',
     6               'Er','Tm','Yb','Lu',
     6'Hf','Ta',' W','Re','Os','Ir','Pt','Au','Hg',
     6          'Tl','Pb','Bi','Po','At','Rn',
     7'Fr','Ra','Ac'/

      inquire(file='GEOMETRY.xyz',exist=lex)
      if(.not.lex)then
       write(6,*)'GEOMETRY.xyz of CPMD not found'
       stop
      endif 
c
      OPEN(8,FILE='GEOMETRY.xyz')
      read(8,*)NAT
      write(6,*)NAT,' atoms'
      if(NAT.gt.NT0)then
       write(6,*)'Too many atoms'
       stop
      endif 
      read(8,*)
      do 1 IA=1,nat
1     read(8,800)AT(IA),(X(IX,IA),IX=1,3)
800   format(1x,a2,3f20.12)
      close(8)
c
      OPEN(15,FILE='FILE.X')
      WRITE(15,*)' FILE.X from CPMD output'
      WRITE(15,*)NAT
      DO 777 IA=1,NAT
      K=0
      do 2 ii=1,89
2     IF(AT(IA).EQ.atsy(ii))K=ii
      IF(K.EQ.0)WRITE(6,*)' Unknown atom ',AT(IA)
777   WRITE(15,1501)K,(X(IX,IA),IX=1,3)
1501  FORMAT(I4,3F15.8,' 0 0 0 0 0 0 0 0.0')
      CLOSE(15)
      WRITE(6,*)' X written into FILE.X'

      inquire(file='HESSIAN',exist=lex)
      if(.not.lex)then
       write(6,*)'HESSIAN of CPMD not found'
       stop
      endif 

      OPEN(8,FILE='HESSIAN')
      read(8,*)N
      read(8,*)((F(IX,IY),IX=1,N),IY=1,N)
      close(8)
      OPEN(20,FILE='FILE.FC')
      CALL WRITEFF(MX3,N,F)
      CLOSE(20)
      WRITE(6,*)' FF written into FILE.FC'

      STOP
      END
 
      SUBROUTINE WRITEFF(MX3,N,FCAR)
      IMPLICIT INTEGER*4 (I-N) 
      IMPLICIT REAL*8(A-H,O-Z)
      DIMENSION FCAR(MX3,N)
C     CONST=4.359828/0.5291772**2
      CONST=1.0d0
      DO 6 I=1,N
      DO 6 J=1,N
6     FCAR(I,J)=FCAR(I,J)/CONST
      N1=1
1     N3=N1+4
      IF(N3.GT.N)N3=N
      DO 130 LN=N1,N
130   WRITE(20,17)LN,(FCAR(LN,J),J=N1,MIN(LN,N3))
      N1=N1+5
      IF(N3.LT.N)GOTO 1
17    FORMAT(I4,5D14.6)
      RETURN
      END