vasp_tools/visualize_contcar

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