'Program to calculate crystal orientation resulting from rotation about ' lab axes, expressed as total rotation from std position in standard ' rotx, roty, rotz. f$="#####.###" 5 pi=4*ATN(1) 6 delta=.001:'closeness to Ewald's sphere within which spot will be observed 7 'f$="###. ###. ###. #.##### #.##### #.##### ####.## ####.##" 'F2$="###.######" 9 distance=100 :'mm xtal to film 10 lambda=1.08:'angstroms 11 lamdainv=1/lambda 12 a(1)=221.4:a(2)=a(1):a(3)=386.95:'angstroms 12 'a(1)=212:a(2)=a(1):a(3)=352:'angstroms 13 gamma(1)=pi/3:gamma(2)=pi/2:gamma(3)=pi/2 14 'for i=1 to 3:b(i)=1/a(i):next i 16 b(1,1)=COS(pi/6)/(a(1)*COS(pi/6)):b(2,1)=SIN(pi/6)/(a(1)*COS(pi/6)):b(3,1)=0 ' 1 = a . a* = |a| |a*|cos 30¡ (a*, b* sep by 60¡, a orthog b*, -> a,a* sep by 30¡) ' |a|=212 ' |a*|=1/(212 cos 30¡) 'take x axis halfway between a*, b*, i.e. both 30¡ from x axis. Component along x is |a*|cos 30¡. 'Component along y is ±|a*|sin(30¡) 'gives 4.72, 2.72, 0 E-3 reciprocal angstroms 17 b(1,2)=b(1,1):b(2,2)=-b(2,1):b(3,2)=0 18 b(1,3)=0:b(2,3)=0:b(3,3)=1/a(3):'0,0, 2.84E-3 rec. angstrom '19 FOR I=1 TO 3: FOR j=1 TO 3:PRINT b(I,j);:NEXT j:PRINT:NEXT I:PRINT 50 'read rot angles: '52 READ rotx,roty,rotz 'PRINT USING f$;rotx, roty,rotz DATA -245.49, 61.32,0 INPUT "standard rotation angles for current position";rotx, roty,rotz 56 rotx=rotx*pi/180:roty=roty*pi/180:rotz=rotz*pi/180 58 rmx(1,1)=1:rmx(2,2)=COS(rotx):rmx(3,2)=SIN(rotx):rmx(2,3)=-rmx(3,2):rmx(3,3)=rmx(2,2) 60 rmy(2,2)=1:rmy(1,1)=COS(roty):rmy(3,1)=SIN(roty):rmy(1,3)=-rmy(3,1):rmy(3,3)=rmy(1,1) 62 rmz(3,3)=1:rmz(1,1)=COS(rotz):rmz(2,1)=SIN(rotz):rmz(1,2)=-rmz(2,1):rmz(2,2)=rmz(1,1) 'multiply the matrices together, put result in prod(). FOR i=1 TO 3:FOR K=1 TO 3:temp(i,K)=0 FOR j=1 TO 3:temp(i,K)=temp(i,K)+rmy(i,j)*rmz(j,K) NEXT j:NEXT K:NEXT i FOR i=1 TO 3:FOR K=1 TO 3:prod(i,K)=0 FOR j=1 TO 3:prod(i,K)=prod(i,K)+rmx(i,j)*temp(j,K) NEXT j:NEXT K:NEXT i postrot: INPUT"further rotate about x,y, or z";a$ INPUT"how many degrees";x x=x*pi/180 IF UCASE$(a$)="X" THEN GOSUB rotx IF UCASE$(a$)="Y" THEN GOSUB roty IF UCASE$(a$)="Z" THEN GOSUB rotz mult: FOR i=1 TO 3:FOR K=1 TO 3:temp(i,K)=0 FOR j=1 TO 3:temp(i,K)=temp(i,K)+rmx(i,j)*prod(j,K) NEXT j:NEXT K:NEXT i '932 FOR i=1 TO 3:FOR j=1 TO 3:PRINT prod(i,j),:NEXT j:PRINT:NEXT i:PRINT FOR i=1 TO 3: FOR j=1 TO 3: prod(i,j)=temp(i,j):NEXT j:NEXT i '926 FOR i=1 TO 3:FOR j=1 TO 3:PRINT rmz(i,j),:NEXT j:PRINT:NEXT i:PRINT '928 FOR i=1 TO 3:FOR j=1 TO 3:PRINT rmy(i,j),:NEXT j:PRINT:NEXT i:PRINT '930 FOR i=1 TO 3:FOR j=1 TO 3:PRINT rmx(i,j),:NEXT j:PRINT:NEXT i:PRINT '932 FOR i=1 TO 3:FOR j=1 TO 3:PRINT prod(i,j),:NEXT j:PRINT:NEXT i:PRINT Factor:' R into 3 rotation matrices x=-prod(1,3):GOSUB arcsin:b=y:'angle of rotation about y IF b>pi THEN b=b-2*pi a=ATN(-prod(2,3)/prod(3,3)) IF prod(3,3)/COS(b)<0 THEN a=a+pi IF a>pi THEN a=a-2*pi c=ATN(-prod(1,2)/prod(1,1)) IF prod(1,1)/COS(b)<0 THEN c=c+pi IF c>pi THEN c=c-2*pi f2$="rotation about x,y,z:####.### ####.### ####.###" PRINT USING f2$; 180*a/pi,180*b/pi,180*c/pi GOTO postrot rotx: rmx(1,2)=0:rmx(1,3)=0:rmx(2,1)=0:rmx(3,1)=0: rmx(1,1)=1 rmx(2,2)=COS(x):rmx(3,2)=SIN(x):rmx(2,3)=-rmx(3,2):rmx(3,3)=rmx(2,2) RETURN roty: rmx(1,2)=0:rmx(3,2)=0:rmx(2,1)=0:rmx(2,3)=0: rmx(2,2)=1 rmx(1,1)=COS(x):rmx(3,1)=SIN(x):rmx(1,3)=-rmx(3,1):rmx(3,3)=rmx(1,1) RETURN rotz: rmx(1,3)=0:rmx(2,3)=0:rmx(3,1)=0:rmx(3,2)=0:rmx(3,3)=1 rmx(1,1)=COS(x):rmx(2,1)=SIN(x):rmx(1,2)=-rmx(2,1):rmx(2,2)=rmx(1,1) RETURN arcsin: IF ABS(x)>1.01 THEN PRINT"invalid arg for arcsin!";x:STOP IF ABS(x)>=1 THEN y=pi/2:IF x<0 THEN y=-y:RETURN y=ATN(x/SQR(1-x*x)) RETURN arccos: IF ABS(x)>1.01 THEN PRINT"invalid arg for arccos!";x:STOP IF ABS(x)>=1 THEN y=pi/2:IF x<0 THEN y=-y:RETURN IF x=0 THEN y=pi/2:RETURN y=ATN(SQR(1-x*x)/x) RETURN notused: 'INPUT"Enter further rotation about x, y, and z (z then y then x)";rotx,roty,rotz x=x*pi/180 'rotx=rotx*pi/180:roty=roty*pi/180:rotz=rotz*pi/180 rmx(1,1)=1:rmx(2,2)=COS(rotx):rmx(3,2)=SIN(rotx):rmx(2,3)=-rmx(3,2):rmx(3,3)=rmx(2,2) rmy(2,2)=1:rmy(1,1)=COS(roty):rmy(3,1)=SIN(roty):rmy(1,3)=-rmy(3,1):rmy(3,3)=rmy(1,1) rmz(3,3)=1:rmz(1,1)=COS(rotz):rmz(2,1)=SIN(rotz):rmz(1,2)=-rmz(2,1):rmz(2,2)=rmz(1,1) FOR i=1 TO 3:FOR K=1 TO 3:temp(i,K)=0 FOR j=1 TO 3:temp(i,K)=temp(i,K)+rmz(i,j)*prod(j,K) NEXT j:NEXT K:NEXT i FOR i=1 TO 3:FOR K=1 TO 3:prod(i,K)=0 FOR j=1 TO 3:prod(i,K)=prod(i,K)+rmy(i,j)*temp(j,K) NEXT j:NEXT K:NEXT i FOR i=1 TO 3:FOR K=1 TO 3:temp(i,K)=0 FOR j=1 TO 3:temp(i,K)=temp(i,K)+rmx(i,j)*prod(j,K) NEXT j:NEXT K:NEXT i FOR i=1 TO 3: FOR j=1 TO 3 prod(i,j)=temp(i,j):NEXT j:NEXT i END 926 FOR i=1 TO 3:FOR j=1 TO 3:PRINT rmz(i,j),:NEXT j:PRINT:NEXT i:PRINT 928 FOR i=1 TO 3:FOR j=1 TO 3:PRINT rmy(i,j),:NEXT j:PRINT:NEXT i:PRINT 930 FOR i=1 TO 3:FOR j=1 TO 3:PRINT rmx(i,j),:NEXT j:PRINT:NEXT i:PRINT 932 FOR i=1 TO 3:FOR j=1 TO 3:PRINT prod(i,j),:NEXT j:PRINT:NEXT i:PRINT OPEN"o",3,"rmat 2" PRINT#3,"rot matrix derived" FOR i=1 TO 3:FOR j=1 TO 3:PRINT#3, USING f2$;prod(i,j),:NEXT:PRINT#3,:NEXT PRINT#3,"RXR:" FOR i=1 TO 3 FOR j=1 TO 3 x=0:FOR K=1 TO 3 x=x+prod(i,K)*prod(j,K) NEXT K PRINT#3, USING f2$; x, NEXT j:PRINT#3, NEXT i CLOSE#3 946 RETURN rotation: 950 FOR i=1 TO 3:xlab(i)=0 952 FOR j=1 TO 3:xlab(i)=xlab(i)+prod(i,j)*x(j) 954 NEXT j:NEXT i 956 FOR i=1 TO 3:x(i)=xlab(i):NEXT i 960 RETURN