-
Notifications
You must be signed in to change notification settings - Fork 0
/
showvib.vmd
131 lines (118 loc) · 4.07 KB
/
showvib.vmd
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
proc drawarrow {atmrange fragdx fragdy fragdz {scl 1} {rad 0.2} {showgoc 1}} {
#This script is used to draw arrow
#Determine arrow center
set sel [atomselect top $atmrange]
set cen [measure center $sel]
set cenx [lindex $cen 0]
set ceny [lindex $cen 1]
set cenz [lindex $cen 2]
if {$showgoc==1} {puts "Geometry center: $cenx $ceny $cenz"}
#Scale vector
set fragdx [expr $fragdx*$scl]
set fragdy [expr $fragdy*$scl]
set fragdz [expr $fragdz*$scl]
#Draw arrow
set body 0.75
set begx [expr $cenx-$fragdx/2]
set begy [expr $ceny-$fragdy/2]
set begz [expr $cenz-$fragdz/2]
set endx [expr $cenx+$fragdx*$body-$fragdx/2]
set endy [expr $ceny+$fragdy*$body-$fragdy/2]
set endz [expr $cenz+$fragdz*$body-$fragdz/2]
draw cylinder "$begx $begy $begz" "$endx $endy $endz" radius $rad filled yes resolution 20
set endx2 [expr $cenx+$fragdx/2]
set endy2 [expr $ceny+$fragdy/2]
set endz2 [expr $cenz+$fragdz/2]
draw cone "$endx $endy $endz" "$endx2 $endy2 $endz2" radius [expr $rad*2.5] resolution 20
}
proc vib {filename imode {sclfac 5} {rad 0.05}} {
#This script is used to draw displacement vectors of vibration generated by Gaussian
#filename.xyz, filename.out should be presented in current folder
color Display Background white
axes location off
display depthcue off
display rendermode GLSL
light 2 off
light 3 on
color Element C silver
color change rgb 6 0.800000 0.800000 0.800000
material change specular Opaque 0.200000
material change transmode EdgyGlass 1.000000
material change specular EdgyGlass 0.200000
proc helpp {} {
puts "No Available Commands"
}
foreach i [molinfo list] {
mol delete $i
}
mol new $filename.xyz
mol modstyle 0 top CPK 0.8 0.3 22.0 22.0
mol modcolor 0 top Element
set ilinear 0
if {[info exist x]} {unset x}
set natm [molinfo top get numatoms]
if {$ilinear==0} {set nmode [expr 3*$natm-6]}
if {$ilinear==1} {set nmode [expr 3*$natm-5]}
puts "Number of modes to load: $nmode"
# Load normal coordinates
set myfile [open $filename.out r]
while { [gets $myfile line] >=0 } {
if {[string first " and normal coordinates:" $line] != -1} {break}
}
set ipos 0
set nframe [expr ceil($nmode/3.0)]
for {set iframe 1} {$iframe<=$nframe} {incr iframe} {
while { [gets $myfile line] >=0 } {
if {[string first "Atom AN" $line] != -1} {break}
}
set nleft [expr $nmode-$ipos]
set m1 [expr $ipos+1]
set m2 [expr $ipos+2]
set m3 [expr $ipos+3]
#puts "Loading frame: $iframe left: $nleft"
if {$nleft>=3} {
for {set iatm 1} {$iatm<=$natm} {incr iatm} {
gets $myfile line
scan $line "%d %d %f %f %f %f %f %f %f %f %f" itmp jtmp \
x($iatm,$m1) y($iatm,$m1) z($iatm,$m1) \
x($iatm,$m2) y($iatm,$m2) z($iatm,$m2) \
x($iatm,$m3) y($iatm,$m3) z($iatm,$m3)
}
incr ipos 3
} elseif {$nleft==2} {
for {set iatm 1} {$iatm<=$natm} {incr iatm} {
gets $myfile line
scan $line "%d %d %f %f %f %f %f %f" itmp jtmp \
x($iatm,$m1) y($iatm,$m1) z($iatm,$m1) \
x($iatm,$m2) y($iatm,$m2) z($iatm,$m2)
}
incr ipos 2
} elseif {$nleft==1} {
for {set iatm 1} {$iatm<=$natm} {incr iatm} {
gets $myfile line
scan $line "%d %d %f %f %f" itmp jtmp x($iatm,$m1) y($iatm,$m1) z($iatm,$m1)
#puts "$x($iatm,$m1) $y($iatm,$m1) $z($iatm,$m1)"
}
incr ipos
}
}
close $myfile
#Print all loaded data
if {0} {
puts "Loaded normal mode:"
for {set imode 1} {$imode<=$nmode} {incr imode} {
puts "Mode $imode:"
for {set iatm 1} {$iatm<=$natm} {incr iatm} {
puts "Atom $iatm: $x($iatm,$imode) $y($iatm,$imode) $z($iatm,$imode)"
}
}
}
#Define command for drawing normal coordinate and print loaded data
draw delete all
draw color yellow
puts "Load Mode $imode"
for {set i 1} {$i<=$natm} {incr i} {
drawarrow "serial $i" $x($i,$imode) $y($i,$imode) $z($i,$imode) $sclfac $rad 0
puts "Atom $iatm: $x($i,$imode) $y($i,$imode) $z($i,$imode)"
}
}