This shows you the differences between two versions of the page.
scripts:modxyz [2009/11/27 12:45] prokop |
scripts:modxyz [2011/02/18 13:13] |
||
---|---|---|---|
Line 1: | Line 1: | ||
- | ==== usage ==== | ||
- | * Copy wf2state to bin directory | ||
- | * copy particular cel.lvs into same directory as your *.xyz file | ||
- | * run with *.xyz file name as 1st parameter and shift by lattice vectors da, db, dc, example: | ||
- | <code> modxyz scan.xyz 0.5 0.5 0.0 </code> | ||
- | * script generate ''scan.xyz.s.xyz'' file | ||
- | |||
- | ==== how it works ==== | ||
- | * compute inverse matrix of lattice | ||
- | * transform geometry by the inverse lattice (make it rectangular) | ||
- | * shift by da,db,dz in this rectangularized space | ||
- | * modulo by cel size a,b,c | ||
- | * transform back by forward lattice matrix | ||
- | |||
- | ==== modxyz: ==== | ||
- | <code> | ||
- | #! /bin/bash | ||
- | |||
- | #input="scan-.xyz" | ||
- | input=$1 | ||
- | output=$input.s.xyz" | ||
- | |||
- | da=$2 | ||
- | db=$3 | ||
- | dc=$4 | ||
- | |||
- | p=`cat $input | head -1` | ||
- | let p=$p+2 | ||
- | |||
- | awk -v p=$p -v da=$da -v db=$db -v dc=$dc ' | ||
- | function deter(a){ | ||
- | det= a[1,1]*(a[3,3]*a[2,2]-a[3,2]*a[2,3])-a[2,1]*(a[3,3]*a[1,2]-a[3,2]*a[1,3])+a[3,1]*(a[2,3]*a[1,2]-a[2,2]*a[1,3]) | ||
- | return det | ||
- | } | ||
- | |||
- | function invert(a,b){ | ||
- | det = deter(a) | ||
- | |||
- | b[1,1] = (a[3,3]*a[2,2]-a[3,2]*a[2,3]) / det | ||
- | b[1,2] = -(a[3,3]*a[1,2]-a[3,2]*a[1,3]) /det | ||
- | b[1,3] = (a[2,3]*a[1,2]-a[2,2]*a[1,3]) /det | ||
- | |||
- | b[2,1] = -(a[3,3]*a[2,1]-a[3,1]*a[2,3]) /det | ||
- | b[2,2] = (a[3,3]*a[1,1]-a[3,1]*a[1,3]) /det | ||
- | b[2,3] = -(a[2,3]*a[1,1]-a[2,1]*a[1,3]) /det | ||
- | |||
- | b[3,1] = (a[3,2]*a[2,1]-a[3,1]*a[2,2]) /det | ||
- | b[3,2] = -(a[3,2]*a[1,1]-a[3,1]*a[1,2]) /det | ||
- | b[3,3] = (a[2,2]*a[1,1]-a[2,1]*a[1,2]) /det | ||
- | } | ||
- | |||
- | function Mvec(M,v,out){ | ||
- | out[1]=v[1]*M[1,1]+v[2]*M[1,2]+v[3]*M[1,3] | ||
- | out[2]=v[1]*M[2,1]+v[2]*M[2,2]+v[3]*M[2,3] | ||
- | out[3]=v[1]*M[3,1]+v[2]*M[3,2]+v[3]*M[3,3] | ||
- | } | ||
- | |||
- | function Mvec2(M,v,out){ | ||
- | out[1]=v[1]*M[1,1]+v[2]*M[2,1]+v[3]*M[3,1] | ||
- | out[2]=v[1]*M[1,2]+v[2]*M[2,2]+v[3]*M[3,2] | ||
- | out[3]=v[1]*M[1,3]+v[2]*M[2,3]+v[3]*M[3,3] | ||
- | } | ||
- | |||
- | |||
- | |||
- | BEGIN{ | ||
- | getline < "cel.lvs" | ||
- | A[1,1]=$1;A[1,2]=$2;A[1,3]=$3 | ||
- | getline < "cel.lvs" | ||
- | A[2,1]=$1;A[2,2]=$2;A[2,3]=$3 | ||
- | getline < "cel.lvs" | ||
- | A[3,1]=$1;A[3,2]=$2;A[3,3]=$3 | ||
- | |||
- | invert(A,B) | ||
- | } | ||
- | |||
- | #(NR==1){print $0} | ||
- | #(NR==2){print $0} | ||
- | |||
- | (NR>0){ if ($4!=""){ | ||
- | v[1]=$2;v[2]=$3;v[3]=$4; | ||
- | Mvec2(B,v,v_) | ||
- | |||
- | v_[1]=(v_[1]+da)%1.0 | ||
- | v_[2]=(v_[2]+db)%1.0 | ||
- | v_[3]=(v_[3]+dc)%1.0 | ||
- | |||
- | #print $1," ",v_[1]," ",v_[2]," ",v_[3] | ||
- | Mvec2(A,v_,v__) | ||
- | print $1," ",v__[1]," ",v__[2]," ",v__[3] | ||
- | }else{print $0} | ||
- | } | ||
- | |||
- | END{ | ||
- | } | ||
- | |||
- | ' $input > $output | ||
- | </code> |