support bonds for more than one atom
[physik/posic.git] / vasp_tools / visualize
1 #!/bin/sh
2
3 #
4 # visualization script
5 # author: frank.zirkelbach@physik.uni-augsburg.de
6 #
7
8 # help function
9 draw_cyl() {
10         cat >> temp.pov <<-EOF
11 cylinder {
12 <$1, $3, $2>, <$4, $6, $5>, 0.05
13 pigment { color White }
14 }
15 EOF
16 }
17 draw_bond() {
18         cat >> temp.pov <<-EOF
19 cylinder {
20 <$1, $3, $2>, <$4, $6, $5>, $7
21 pigment { color Blue }
22 }
23 EOF
24 }
25
26 # defaults
27
28 lc=5.429
29 directory="doesnt_exist____for_sure"
30 width="640"
31 height="480"
32 radius="0.6"
33 x0="-0.6"; y0="-0.6"; z0="-0.6";
34 x1="0.6"; y1="0.6"; z1="0.6";
35 cx=""; cy=""; cz="";
36 lx="0"; ly="-100"; lz="100";
37 ortographic=""
38 bx0=""; by0=""; bz0="";
39 bx1=""; by1=""; bz1="";
40 bcr="0.1";
41 clx="0"; cly="0"; clz="0";
42 extra=0
43 displace=""
44 mirror=0
45 mx1=0; mx2=0; mx3=0;
46 my1=0; my2=0; my3=0;
47 mz1=0; mz2=0; mz3=0;
48 ab=0
49
50 # parse argv
51
52 while [ "$1" ]; do
53         case "$1" in
54                 -d)             directory=$2;           shift 2;;
55                 -w)             width=$2;               shift 2;;
56                 -h)             height=$2;              shift 2;;
57                 -r)             radius=$2;              shift 2;;
58                 -nll)           x0=$2; y0=$3; z0=$4;    shift 4;;
59                 -fur)           x1=$2; y1=$3; z1=$4;    shift 4;;
60                 -c)             cx=$2; cy=$3; cz=$4;    shift 4;;
61                 -L)             clx=$2; cly=$3; clz=$4; shift 4;;
62                 -l)             lx=$2; ly=$3; lz=$4;    shift 4;;
63                 -o)             ortographic=1;          shift 1;;
64                 -b)             bx0=$2; by0=$3; bz0=$4;
65                                 bx1=$5; by1=$6; bz1=$7; shift 7;;
66                 -B)             bcr=$2;                 shift 2;;
67                 -C)             lc=$2;                  shift 2;;
68                 -e)             extra=1;                shift 1;;
69                 -D)             displace=$2;            shift 2;;
70                 -m)             mx1=$2; mx2=$3; mx3=$4;
71                                 my1=$5; my2=$6; my3=$7;
72                                 mz1=$8; mz2=$9; mz3=${10};
73                                 mirror=1;               shift 10;;
74                 -A)             ab=$2;                  shift 2;
75                                 ((cnt=1))
76                                 while [ $cnt -le $ab ]; do
77                                         anr[$cnt]=$1;   shift 1;
78                                         ((cnt+=1))
79                                 done
80                                 cutoff=$1;              shift 1;;
81                 *)
82                                 echo "options:"
83                                 echo "########"
84                                 echo "directory to progress:"
85                                 echo "  -d <directory> (mandatory)"
86                                 echo "png dim:"
87                                 echo "  -w <width>"
88                                 echo "  -h <height>"
89                                 echo "atom size:"
90                                 echo "  -r <radius>"
91                                 echo "  -B <bond cylinder radius>"
92                                 echo "unit cell:"
93                                 echo "  -C <lattice constant>"
94                                 echo "  -m <dimX> <dimY> <dimZ> (mirror atoms)"
95                                 echo "visualization volume:"
96                                 echo "  -nll <x> <y> <z> (near lower left)"
97                                 echo "  -fur <x> <y> <z> (far upper right)"
98                                 echo "  -o (ortographic)"
99                                 echo "bounding box:"
100                                 echo "  -b <x0> <y0> <z0> <x1> <y1> <z1>"
101                                 echo "povray:"
102                                 echo "  -c <x> <y> <z> (camera position)"
103                                 echo "  -L <x> <y> <z> (camera look)"
104                                 echo "  -l <x> <y> <z> (light source)"
105                                 echo "bonds:"
106                                 echo "  -ab <atom number> <cutoff> (auto bonds)"
107                                 exit 1;;
108         esac
109 done
110
111 # calculation from lattic eunits to angstroms
112
113 [ "$lc" = "sic" ] && lc=4.359
114 [ "$lc" = "si" ] && lc=5.480
115 [ "$lc" = "c" ] && lc=3.566
116
117 #offset=`echo 0.125 \* $lc | bc`
118 offset=0.0
119
120 x0=`echo $x0 \* $lc + $offset | bc`
121 y0=`echo $y0 \* $lc + $offset | bc`
122 z0=`echo $z0 \* $lc + $offset | bc`
123 x1=`echo $x1 \* $lc + $offset | bc`
124 y1=`echo $y1 \* $lc + $offset | bc`
125 z1=`echo $z1 \* $lc + $offset | bc`
126
127 mx1=`echo $mx1 \* $lc + $offset | bc`
128 my1=`echo $my1 \* $lc + $offset | bc`
129 mz1=`echo $mz1 \* $lc + $offset | bc`
130 mx2=`echo $mx2 \* $lc + $offset | bc`
131 my2=`echo $my2 \* $lc + $offset | bc`
132 mz2=`echo $mz2 \* $lc + $offset | bc`
133 mx3=`echo $mx3 \* $lc + $offset | bc`
134 my3=`echo $my3 \* $lc + $offset | bc`
135 mz3=`echo $mz3 \* $lc + $offset | bc`
136
137 clx=`echo $clx \* $lc + $offset | bc`
138 cly=`echo $cly \* $lc + $offset | bc`
139 clz=`echo $clz \* $lc + $offset | bc`
140
141 if [ -n "$cx" -a -n "$cy" -a -n "$cz" ]; then
142         cx=`echo $cx \* $lc + $offset | bc`
143         cy=`echo $cy \* $lc + $offset | bc`
144         cz=`echo $cz \* $lc + $offset | bc`
145 fi
146
147 if [ -n "$bx0" ]; then
148         bx0=`echo $bx0 \* $lc + $offset | bc`
149         by0=`echo $by0 \* $lc + $offset | bc`
150         bz0=`echo $bz0 \* $lc + $offset | bc`
151         bx1=`echo $bx1 \* $lc + $offset | bc`
152         by1=`echo $by1 \* $lc + $offset | bc`
153         bz1=`echo $bz1 \* $lc + $offset | bc`
154 fi
155
156 # povray command
157
158 POVRAY="povray -W${width} -H${height} -d" 
159
160 # convert options
161
162 COPTS="-font /usr/share/fonts/truetype/ttf-dejavu/DejaVuSans.ttf"
163 COPTS="$COPTS -depth 8 -fill white -stroke blue -pointsize 24"
164
165 # do it ...
166
167 if [ -d $directory ]; then
168         filesource=$directory/atomic_conf_*.xyz
169 fi
170
171 if [ -f $directory ]; then
172         filesource=$directory
173 fi
174
175 for file in $filesource; do
176
177         cat > temp.pov <<-EOF
178 #include "colors.inc"
179 #include "textures.inc"
180 #include "shapes.inc"
181 #include "glass.inc"
182 #include "metals.inc"
183 #include "woods.inc"
184 #include "stones.inc"
185 EOF
186
187         # meta info
188         count=`grep '# \[P\]' $file | awk '{ print $3 }'`
189         time=`grep '# \[P\]' $file | awk '{ print $4 }'`
190         camloc=`grep '# \[P\]' $file | awk '{ print $5 }'`
191         [ -n "$cx" -a -n "$cy" -a -n "$cz" ] && camloc="<$cx,$cz,$cy>"
192
193         # atoms
194 #echo "-----> $x0 $y0 $z0 $x1 $y1 $z1 <-----"
195         if [ -n "$x0" ]; then
196                 export x0 y0 z0 x1 y1 z1 radius extra
197                 export mx1 mx2 mx3
198                 export my1 my2 my3
199                 export mz1 mz2 mz3
200                 export mirror
201                 cat $file | grep -v '#' | awk '\
202                 BEGIN {
203                         x0=ENVIRON["x0"]; y0=ENVIRON["y0"]; z0=ENVIRON["z0"];
204                         x1=ENVIRON["x1"]; y1=ENVIRON["y1"]; z1=ENVIRON["z1"];
205                         mx1=ENVIRON["mx1"]; my1=ENVIRON["my1"];
206                         mz1=ENVIRON["mz1"]; mx2=ENVIRON["mx2"];
207                         my2=ENVIRON["my2"]; mz2=ENVIRON["mz2"];
208                         mx3=ENVIRON["mx3"]; my3=ENVIRON["my3"];
209                         mz3=ENVIRON["mz3"]; mirror=ENVIRON["mirror"];
210                         radius=ENVIRON["radius"]; extra=ENVIRON["extra"];
211                         nx=0; ny=0; nz=0;
212                 }
213                 {
214                         if(($2>=x0)&&($3>=y0)&&($4>=z0)&&\
215                            ($2<=x1)&&($3<=y1)&&($4<=z1)) {
216                                 print "sphere { <"$2","$4","$3">, "radius" ";
217                                 if(extra)
218                 print "texture { pigment { color rgb<"$6/4.4",0,"1-$6/4.4"> } ";
219                                 else
220                 print "texture { pigment { color "$5" } ";
221                                 print "finish { phong 1 metallic } } }";
222                         }
223
224                         if(mirror) {
225
226                         nx=$2-mx1; ny=$3-mx2; nz=$4-mx3;
227                         if((nx>=x0)&&(ny>=y0)&&(nz>=z0)&&\
228                            (nx<=x1)&&(ny<=y1)&&(nz<=z1)) {
229 print "// translation: "$2/5.480" "$3/5.480" "$4/5.480" -> "nx/5.480" "ny/5.480" "nz/5.480" ...";
230                                 print "sphere { <"nx","nz","ny">, "radius" ";
231                                 if(extra)
232                 print "texture { pigment { color rgb<"$6/4.4",0,"1-$6/4.4"> } ";
233                                 else
234                 print "texture { pigment { color "$5" } ";
235                                 print "finish { phong 1 metallic } } }";
236                         }
237                         nx=$2-my1; ny=$3-my2; nz=$4-my3;
238                         if((nx>=x0)&&(ny>=y0)&&(nz>=z0)&&\
239                            (nx<=x1)&&(ny<=y1)&&(nz<=z1)) {
240                                 print "sphere { <"nx","nz","ny">, "radius" ";
241                                 if(extra)
242                 print "texture { pigment { color rgb<"$6/4.4",0,"1-$6/4.4"> } ";
243                                 else
244                 print "texture { pigment { color "$5" } ";
245                                 print "finish { phong 1 metallic } } }";
246                         }
247                         nx=$2-mz1; ny=$3-mz2; nz=$4-mz3;
248                         if((nx>=x0)&&(ny>=y0)&&(nz>=z0)&&\
249                            (nx<=x1)&&(ny<=y1)&&(nz<=z1)) {
250                                 print "sphere { <"nx","nz","ny">, "radius" ";
251                                 if(extra)
252                 print "texture { pigment { color rgb<"$6/4.4",0,"1-$6/4.4"> } ";
253                                 else
254                 print "texture { pigment { color "$5" } ";
255                                 print "finish { phong 1 metallic } } }";
256                         }
257
258                         }
259
260                 }' >> temp.pov
261         else
262                 cat $file | grep -v '#' | while read name x y z color temp; do
263                         cat >> temp.pov <<-EOF
264 sphere {
265 <$x, $z, $y>, $radius
266 texture {
267 pigment { color $color }
268 finish {
269 phong 1
270 metallic
271 }
272 }
273 }
274 EOF
275                 done
276         fi
277
278         # boundaries
279         if [ -z "$bx0" ]; then
280
281         #if [ -z "$x0" ]; then
282
283         cat $file | grep '# \[D\]' | while read foo bar x1 y1 z1 x2 y2 z2 ; do
284                 draw_cyl $x1 $y1 $z1 $x2 $y2 $z2 0.05
285         done
286
287         #else
288
289                 # manually drawing the 3x4 boundaries ...
290 #               draw_cyl $x0 $y0 $z0 $x1 $y0 $z0
291 #               draw_cyl $x0 $y0 $z0 $x0 $y1 $z0
292 #               draw_cyl $x1 $y1 $z0 $x1 $y0 $z0
293 #               draw_cyl $x0 $y1 $z0 $x1 $y1 $z0
294
295 #               draw_cyl $x0 $y0 $z1 $x1 $y0 $z1
296 #               draw_cyl $x0 $y0 $z1 $x0 $y1 $z1
297 #               draw_cyl $x1 $y1 $z1 $x1 $y0 $z1
298 #               draw_cyl $x0 $y1 $z1 $x1 $y1 $z1
299
300 #               draw_cyl $x0 $y0 $z1 $x0 $y0 $z0
301 #               draw_cyl $x0 $y1 $z1 $x0 $y1 $z0
302 #               draw_cyl $x1 $y0 $z1 $x1 $y0 $z0
303 #               draw_cyl $x1 $y1 $z1 $x1 $y1 $z0
304 #       fi
305
306         else
307
308                 # manually drawing the 3x4 boundaries specified by argv ...
309                 draw_cyl $bx0 $by0 $bz0 $bx1 $by0 $bz0
310                 draw_cyl $bx0 $by0 $bz0 $bx0 $by1 $bz0
311                 draw_cyl $bx1 $by1 $bz0 $bx1 $by0 $bz0
312                 draw_cyl $bx0 $by1 $bz0 $bx1 $by1 $bz0
313
314                 draw_cyl $bx0 $by0 $bz1 $bx1 $by0 $bz1
315                 draw_cyl $bx0 $by0 $bz1 $bx0 $by1 $bz1
316                 draw_cyl $bx1 $by1 $bz1 $bx1 $by0 $bz1
317                 draw_cyl $bx0 $by1 $bz1 $bx1 $by1 $bz1
318
319                 draw_cyl $bx0 $by0 $bz1 $bx0 $by0 $bz0
320                 draw_cyl $bx0 $by1 $bz1 $bx0 $by1 $bz0
321                 draw_cyl $bx1 $by0 $bz1 $bx1 $by0 $bz0
322                 draw_cyl $bx1 $by1 $bz1 $bx1 $by1 $bz0
323
324         fi      
325
326         # bonds
327         if [ -n "$bcr" ]; then
328
329                 if [ -z "$x0" ]; then
330
331         cat $file | grep '# \[B\]' | while read foo bar x1 y1 z1 x2 y2 z2 ; do
332                 draw_bond $x1 $z1 $y1 $x2 $z2 $y2 $bcr
333         done
334
335                 else
336
337                 export x0 y0 z0 x1 y1 z1 bcr
338                 cat $file | grep '# \[B\]' | awk '\
339                 BEGIN {
340                         x0=ENVIRON["x0"]; y0=ENVIRON["y0"]; z0=ENVIRON["z0"];
341                         x1=ENVIRON["x1"]; y1=ENVIRON["y1"]; z1=ENVIRON["z1"];
342                         bcr=ENVIRON["bcr"];
343                 }
344                 {
345                         if(($3>=x0)&&($4>=y0)&&($5>=z0)&&\
346                            ($3<=x1)&&($4<=y1)&&($5<=z1)) {
347                                 print "cylinder {";
348                                 print "<"$3","$5","$4">,";
349                                 print "<"$6","$8","$7">, "bcr;
350                                 print "pigment { color Blue }";
351                                 print "}";
352                         }
353                 }' >> temp.pov
354
355                 fi
356         fi
357
358         # auto bonds
359         if [ "$ab" != "0" ]; then
360                 ((cnt=1))
361                 while [ $cnt -le $ab ]; do
362                         anr=${anr[$cnt]}
363
364                 ((tmp=anr+1))
365                 line="`sed -n ${tmp}p $file`"
366                 aX=`echo $line | awk '{ print $2 }'`
367                 aY=`echo $line | awk '{ print $3 }'`
368                 aZ=`echo $line | awk '{ print $4 }'`
369                 export aX aY aZ cutoff mirror
370                 cat $file | grep -v '#' | awk '\
371                 BEGIN {
372                         x=ENVIRON["aX"]; y=ENVIRON["aY"]; z=ENVIRON["aZ"];
373                         co=ENVIRON["cutoff"]; dist=0; bcr=ENVIRON["bcr"];
374                         mx1=ENVIRON["mx1"]; my1=ENVIRON["my1"];
375                         mz1=ENVIRON["mz1"]; mx2=ENVIRON["mx2"];
376                         my2=ENVIRON["my2"]; mz2=ENVIRON["mz2"];
377                         mx3=ENVIRON["mx3"]; my3=ENVIRON["my3"];
378                         mz3=ENVIRON["mz3"]; mirror=ENVIRON["mirror"];
379                 }
380                 {
381                         dist=sqrt((x-$2)^2+(y-$3)^2+(z-$4)^2);
382                         if((dist<=co)&&(dist>0.01)) {
383                                 print "cylinder {";
384                                 print "<"x","z","y">,";
385                                 print "<"$2","$4","$3">, "bcr;
386                                 print "pigment { color Blue }";
387                                 print "}";
388                         }
389
390                         if(mirror) {
391                                 nx=$2-mx1; ny=$3-mx2; nz=$4-mx3;
392                                 dist=sqrt((x-nx)^2+(y-ny)^2+(z-nz)^2)
393                                 if((dist<=co)&&(dist>0.01)) {
394                                         print "cylinder {";
395                                         print "<"x","z","y">,";
396                                         print "<"nx","nz","ny">, "bcr;
397                                         print "pigment { color Blue }";
398                                         print "}";
399                                 }
400                                 nx=$2-my1; ny=$3-my2; nz=$4-my3;
401                                 dist=sqrt((x-nx)^2+(y-ny)^2+(z-nz)^2)
402                                 if((dist<=co)&&(dist>0.01)) {
403                                         print "cylinder {";
404                                         print "<"x","z","y">,";
405                                         print "<"nx","nz","ny">, "bcr;
406                                         print "pigment { color Blue }";
407                                         print "}";
408                                 }
409                                 nx=$2-mz1; ny=$3-mz2; nz=$4-mz3;
410                                 dist=sqrt((x-nx)^2+(y-ny)^2+(z-nz)^2)
411                                 if((dist<=co)&&(dist>0.01)) {
412                                         print "cylinder {";
413                                         print "<"x","z","y">,";
414                                         print "<"nx","nz","ny">, "bcr;
415                                         print "pigment { color Blue }";
416                                         print "}";
417                                 }
418                         }
419
420                 }' >> temp.pov
421
422                         ((cnt+=1))
423                 done
424         fi
425
426         # add camera and light source
427         cat >> temp.pov <<-EOF
428 camera {
429 EOF
430         if [ -n "$ortographic" ]; then  cat >> temp.pov <<-EOF
431 orthographic
432 EOF
433         fi
434         cat >> temp.pov <<-EOF
435 location $camloc
436 look_at <$clx,$clz,$clz>
437 }
438 light_source { <0,100,-100> color White shadowless }
439 EOF
440
441         # mv png
442         $POVRAY temp.pov > /dev/null 2>&1
443         time=`basename $file | awk -F. '{ print $1 }' | awk -F_ '{ print $3 }'`
444         if [ ! -z $displace ]; then
445         convert $COPTS -draw "text 5,20 'd = $displace %'" temp.png temp.png
446         else
447         convert $COPTS -draw "text 5,20 't = $time fs'" temp.png temp.png
448         fi
449         mv temp.png `echo $file | sed 's/\.xyz/\.png/'`
450
451 done
452
453 echo "done"