program gahtables IMPLICIT none real*8 CM,adum character*80 fn character*120 s80,s120 logical l2,lram,lroa,labs,lvcd,lds,lspec,lrs integer*4 i,M,ifa,ifb,io,nb,ic,id,iof, 1nt,nlinepar,iargc real*8,allocatable:: e(:),wh(:),ds(:),rs(:), 1ram180(:),roa180(:),ramdcpi(:),ramdcpii(:),dsh(:), 1roadcpi(:),roadcpii(:) c write(6,960) 960 format(' Extracts anharmonic spectral intensities',/, 1 ' from Gaussian output',/,/, 1 ' Input: G.OUT (or specify filename inline)',/,/, 1 ' Output: DH.TAB harmonic absorption',/, 1 ' DAH.TAB anharmonic absorption',/, 1 ' ROAH.TAB harmonic Raman and ROA',/, 1 ' ROAAH.TAB anharmonic Raman and ROA',/,/) CM=219470.0d0 l2=.false. labs=.false. lvcd=.false. lds=.false. lrs=.false. lram=.false. lroa=.false. lspec=.false. M=0 ifa=0 ifb=0 io=0 nb=0 ic=0 id=0 nt=0 iof=0 nlinepar=iargc() if(nlinepar.gt.0)then call getarg(1,fn) else fn='G.OUT' endif 4 open(9,file='G.OUT',status='old') 1 read(9,900,end=999,err=999)s80 900 format(a120) if(s80(9:49).eq.'QUADRATIC FORCE CONSTANTS IN NORMAL MODES')then call s8 if(.not.l2)then l2=.true. 3 read(9,*,end=88,err=88)adum,adum,adum M=M+1 goto 3 88 close(9) write(6,*)M,' modes ' if(M.gt.0)then allocate(e(M)) do 5 i=1,M 5 e(i)=0.0d0 goto 4 endif call report('No modes found') else do 7 i=1,M read(9,*)e(i),e(i),e(i) 7 e(i)=e(i)/CM endif endif c c G16 output order: c c Anharmonic Infrared Spectroscopy c c 1 Integrated intensity (I) in km.mol^-1 c 2 Dipole strengths (DS) in 10^-40 esu^2.cm^2 c c Anharmonic Raman Spectroscopy (Dynamic) c c 3 SCP(90)z ## Raman activity (RA) in Ang^6 c 4 SCP(180) ## Raman activity (RA) in Ang^6 c 5 DCP_I(180) ## Raman activity (RA) in Ang^6 c c Raman activity (RA) in Ang^6 c c 6 SCP(90)z ## ROA Intensity (ROA) in 10^-4 Ang^6 c 7 SCP(180) ## ROA Intensity (ROA) in 10^-4 Ang^6 c 8 DCP_I(180) ## ROA Intensity (ROA) in 10^-4 Ang^6 c 9 SCP(90)z ##Dimensionless circular intensity difference (CID) c 10 SCP(180) ##Dimensionless circular intensity difference (CID) c 11 DCP_I(180) ##Dimensionless circular intensity difference (CID) if(s80(15:46).eq.'Anharmonic Infrared Spectroscopy')then write(6,900)s80 call af(labs,lvcd,lram,lroa) labs=.true. endif if(s80(17:43).eq.'Anharmonic VCD Spectroscopy')then write(6,900)s80 call af(labs,lvcd,lram,lroa) lvcd=.true. nt=0 endif if(s80(12:50).eq.'Anharmonic Raman Spectroscopy (Dynamic)'.or. 1 s80(11:49).eq.'Anharmonic Raman Spectroscopy (Dynamic)')then write(6,900)s80 call af(labs,lvcd,lram,lroa) lram=.true. id=0 nt=0 endif if(s80(14:46).eq.'Anharmonic Raman Optical Activity')then write(6,900)s80 call af(labs,lvcd,lram,lroa) lroa=.true. id=0 nb=0 endif if(labs.and.s80(9:32).eq.'Dipole strengths (DS) in')then write(6,900)s80 lds=.true. open(11,file='ABS.TXT.SCR') endif if(lvcd.and.s80(9:36).eq.'Rotational strengths (RS) in')then write(6,900)s80 lrs=.true. open(11,file='VCD.TXT.SCR') endif if(s80(17:37).eq.'INFORMATION: SCP(90)z')then id=0 nt=0 if(lroa.or.lram)then write(6,900)s80 if(lram)open(11,file='RAM.SCP90z.TXT.SCR') if(lroa)open(11,file='ROA.SCP90z.TXT.SCR') endif endif if(s80(17:37).eq.'INFORMATION: SCP(180)')then id=0 nt=0 if(lroa.or.lram)then write(6,900)s80 close(11) if(lram)open(11,file='RAM.SCP180.TXT.SCR') if(lroa)open(11,file='ROA.SCP180.TXT.SCR') endif endif if(s80(17:39).eq.'INFORMATION: DCP_I(180)')then id=1 nt=0 if(lroa.or.lram)then write(6,900)s80 close(11) if(lroa)open(11,file='ROA.DCPI.TXT.SCR') if(lram)open(11,file='RAM.DCPI.TXT.SCR') endif endif lspec=lroa.or.lram.or.lds.or.lrs if(lspec)then IF(s80(2:18).eq.'Fundamental Bands'.or. 1 s80(2:10).eq.'Overtones'.or. 2 s80(2:18).eq.'Combination Bands'.or. 2 s80(2:19).eq.'Vibrational States')THEN c ("Vibrational States" for development version) read(9,*) read(9,900)s120 if(s80(2:19).eq.'Vibrational States')THEN iof=-1 else call dof(s120,iof) endif nb=0 if(s80(2:18).eq.'Fundamental Bands')ifb=1 if(s80(2:10).eq.'Overtones')io=1 if(s80(2:18).eq.'Combination Bands')ic=1 if(s80(2:19).eq.'Vibrational States')ifa=1 goto 1 ENDIF IF(ifa.eq.1.or.ifb.eq.1.or.io.eq.1.or.ic.eq.1)THEN if((s80(10:13).eq.') ').or. 1 (ifa.eq.1.and.s80(8:8).eq.' '.and.s80(13:13).ne.' '))then write(11,900)s80 nb=nb+1 nt=nt+1 else if(ifa.eq.1)then write(6,*)nb,' Vibrational States' endif if(ifb.eq.1)then write(6,*)nb,' fundamental modes' ifb=0 endif if(io.eq.1)then write(6,*)nb,' overtones' io=0 endif if(ic.eq.1.or.ifa.eq.1)then write(6,*)nb,' combination bands, ',nt,' total' ic=0 ifa=0 if(lds)then close(11) allocate(wh(nt),ds(nt),dsh(nt)) call suck('ABS.TXT.SCR',nt,25,36,wh,iof,2) call suck('ABS.TXT.SCR',nt,48,63,dsh,iof,3) do 201 i=1,nt 201 dsh(i)=dsh(i)*1.0d-4 call dt('DH.TAB',nt,wh,dsh,dsh) call suck('ABS.TXT.SCR',nt,38,47,wh,iof,2) call suck('ABS.TXT.SCR',nt,64,79,ds,iof,3) do 202 i=1,nt 202 ds(i)=ds(i)*1.0d-4 call dt('DAH.TAB',nt,wh,ds,ds) lds=.false. labs=.false. call system('rm ABS.TXT.SCR') deallocate(wh) endif if(lrs)then close(11) allocate(wh(nt),rs(nt)) call suck('VCD.TXT.SCR',nt,25,36,wh,iof,2) call suck('VCD.TXT.SCR',nt,48,63,rs,iof,3) do 205 i=1,nt 205 rs(i)=rs(i)*1.0d-8 call dt('DH.TAB',nt,wh,dsh,rs) call suck('VCD.TXT.SCR',nt,38,47,wh,iof,2) call suck('VCD.TXT.SCR',nt,64,79,rs,iof,3) do 206 i=1,nt 206 rs(i)=rs(i)*1.0d-8 call dt('DAH.TAB',nt,wh,ds,rs) lrs=.false. lvcd=.false. call system('rm VCD.TXT.SCR') deallocate(wh) endif if(id.eq.1)then if(lram)lram=.false. if(lroa.and.id.eq.1)then lroa=.false. close(11) allocate(wh(nt),ram180(nt),roa180(nt), 1 ramdcpi(nt),roadcpi(nt),ramdcpii(nt),roadcpii(nt)) call suck('RAM.SCP180.TXT.SCR',nt,25,36,wh,iof,2) call suck('RAM.SCP180.TXT.SCR',nt,48,63,ram180,iof,3) call suck('ROA.SCP180.TXT.SCR',nt,48,63,roa180,iof,3) call mw(nt,roa180,ram180,wh) call dr('ROAH.TAB',nt,wh,ram180,roa180) call suck('RAM.DCPI.TXT.SCR',nt,25,36,wh,iof,2) call suck('RAM.DCPI.TXT.SCR',nt,48,63,ramdcpi,iof,3) call suck('ROA.DCPI.TXT.SCR',nt,48,63,roadcpi,iof,3) call mw(nt,roadcpi,ramdcpi,wh) call dr('ROAH_DCPI.TAB',nt,wh,ramdcpi,roadcpi) do 207 i=1,nt roadcpii(i)=0.0d0 207 ramdcpii(i)=ram180(i)-ramdcpi(i) call dr('ROAH_DCPII.TAB',nt,wh,ramdcpii,roadcpii) call suck('RAM.SCP180.TXT.SCR',nt,38,47,wh,iof,2) call suck('RAM.SCP180.TXT.SCR',nt,64,79,ram180,iof,3) call suck('ROA.SCP180.TXT.SCR',nt,64,79,roa180,iof,3) call mw(nt,roa180,ram180,wh) call dr('ROAAH.TAB',nt,wh,ram180,roa180) call suck('RAM.DCPI.TXT.SCR',nt,38,47,wh,iof,2) call suck('RAM.DCPI.TXT.SCR',nt,64,79,ramdcpi,iof,3) call suck('ROA.DCPI.TXT.SCR',nt,64,79,roadcpi,iof,3) call mw(nt,roadcpi,ramdcpi,wh) call dr('ROAAH_DCPI.TAB',nt,wh,ramdcpi,roadcpi) do 208 i=1,nt roadcpii(i)=0.0d0 208 ramdcpii(i)=ram180(i)-ramdcpi(i) call dr('ROAAH_DCPII.TAB',nt,wh,ramdcpii,roadcpii) call system('rm RAM.SCP180.TXT.SCR') call system('rm ROA.SCP180.TXT.SCR') call system('rm RAM.DCPI.TXT.SCR') call system('rm ROA.DCPI.TXT.SCR') call system('rm RAM.SCP90z.TXT.SCR') call system('rm ROA.SCP90z.TXT.SCR') endif endif endif endif ENDIF endif goto 1 999 close(9) END subroutine s8 do 2 i=1,8 2 read(9,*) return end subroutine report(s) character*(*) s write(6,*)s stop end subroutine suck(s,N,i,j,x,iof,ios) implicit none character*(*) s integer*4 iof,N,i,j,k,l,ic,ir,jr,ios real*8 x(*) character*120 st ir=i+iof jr=j+iof open(11,file=s) do 1 k=1,N read(11,100)st 100 format(a120) if(iof.ge.0)then ic=0 do 2 l=ir,jr 2 if(st(l:l).ne.' ')ic=ic+1 if(ic.eq.0)st(jr:jr)='0' read(st(ir:jr),*,end=99,err=99)x(k) else if(ios.eq.2)then read(st,*,end=99,err=99)x(k),x(k) else read(st,*,end=99,err=99)x(k),x(k),x(k) endif endif 1 continue close(11) return 99 write(6,100)st write(6,200)ir,jr,N,k,iof 200 format('Cannot get number between ',i3,' and ',i3,/, 1' N = ',i5,', current: ',i5,' iof: ',i3) stop end subroutine dt(s,N,w,d,r) implicit none character*(*) s integer*4 N,I real*8 w(*),d(*),r(*),pi,po pi=4.0d0*atan(1.0d0) po=dsqrt(2.0d0*pi) open(18,file=s) WRITE (18,501) 501 FORMAT (/,' FREQ. ', 1' D [D^2] ',' R [D^-2]DO R - CO I[km/mol]', 2/,'-----------------------------------------------------------') DO 412 I=1,N 412 WRITE (18,523) I,w(I),d(I),r(I),r(I),po*w(I)*d(I) 523 FORMAT(I5,f16.3,4G15.6) WRITE(18,1819) 1819 FORMAT('-------------------------------------------------') CLOSE(18) write(6,*)s//' written' return end subroutine dr(s,N,w,d,r) implicit none character*(*) s integer*4 N,I real*8 w(*),d(*),r(*),YDX,YDY,CID2,CID1,P1,DX,roa1,doc,roa3, 1ram open(4,file=s) WRITE(4,3304) 3304 FORMAT( 1'MODE FREQ Ramx (180) Ramy Ramtot CID90', 2' CID-X CID180 DEP ROA180 RamanDCPI',/, 3' cm-1 A^4/AMU x10^4', 4' x10^4 x10^4 A^5/AMUx1e4 au.10^4', 5/80(1H-)) do 1 I=1,N YDX=d(I) YDY=d(I) ram=d(I) CID2=0.0d0 CID1=0.0d0 P1=1.0d0 DX=0.0d0 roa1=r(I) doc=r(I) roa3=r(I) 1 WRITE(4,3000)I,w(I),YDX,YDY,ram,CID2,DX,CID1,P1,roa1,doc,roa3 3000 FORMAT(I5,f9.2,6g12.3,f6.3,4g12.4) WRITE(4,1819) 1819 FORMAT('-------------------------------------------------') CLOSE(4) write(6,*)s//' written' RETURN END subroutine af(labs,lvcd,lram,lroa) logical labs,lvcd,lram,lroa labs=.false. lvcd=.false. lram=.false. lroa=.false. return end subroutine dof(s,iof) implicit none character*(*) s integer*4 i,iof do 1 i=1,len(s) if(s(i:i+6).eq.'E(harm)')then iof=i-29 return endif 1 continue call report('cannot see offset') end subroutine mw(nt,roa180,ram180,wh) implicit none integer*4 nt,i real*8 roa180(*),ram180(*),wh(*) do 203 i=1,nt roa180(i)=roa180(i)*wh(i)*0.05932d0 203 ram180(i)=ram180(i)*wh(i)*0.05932d0 return end