34
35
36
37#include "implicit_f.inc"
38
39
40
41#include "mvsiz_p.inc"
42
43
44
45
47 . x1(*), x2(*), x3(*), x4(*), x5(*), x6(*), x7(*), x8(*),
48 . y1(*), y2(*), y3(*), y4(*), y5(*), y6(*), y7(*), y8(*),
49 . z1(*), z2(*), z3(*), z4(*), z5(*), z6(*), z7(*), z8(*),
50 . rx(*) ,ry(*) ,rz(*) ,sx(*) ,sy(*) ,sz(*) ,tx(*) ,ty(*) ,tz(*),
51 . e1x(*),e1y(*),e1z(*),e2x(*),e2y(*),e2z(*),e3x(*),e3y(*),e3z(*)
52
53
54
55#include "vect01_c.inc"
56
57
58
59 INTEGER I
60
62 . x17(mvsiz) , x28(mvsiz) , x35(mvsiz) , x46(mvsiz),
63 . y17(mvsiz) , y28(mvsiz) , y35(mvsiz) , y46(mvsiz),
64 . z17(mvsiz) , z28(mvsiz) , z35(mvsiz) , z46(mvsiz),
65 . a17(mvsiz) , a28(mvsiz) ,
66 . b17(mvsiz) , b28(mvsiz) ,
67 . c17(mvsiz) , c28(mvsiz) ,
68 . det,c1,c2
69
70 DO 10 i=lft,llt
71 x17(i)=x7(i)-x1(i)
72 x28(i)=x8(i)-x2(i)
73 x35(i)=x5(i)-x3(i)
74 x46(i)=x6(i)-x4(i)
75 y17(i)=y7(i)-y1(i)
76 y28(i)=y8(i)-y2(i)
77 y35(i)=y5(i)-y3(i)
78 y46(i)=y6(i)-y4(i)
79 z17(i)=z7(i)-z1(i)
80 z28(i)=z8(i)-z2(i)
81 z35(i)=z5(i)-z3(i)
82 z46(i)=z6(i)-z4(i)
83 10 CONTINUE
84 DO 20 i=lft,llt
85 rx(i)=x17(i)+x28(i)-x35(i)-x46(i)
86 ry(i)=y17(i)+y28(i)-y35(i)-y46(i)
87 rz(i)=z17(i)+z28(i)-z35(i)-z46(i)
88 a17(i)=x17(i)+x46(i)
89 a28(i)=x28(i)+x35(i)
90 b17(i)=y17(i)+y46(i)
91 b28(i)=y28(i)+y35(i)
92 c17(i)=z17(i)+z46(i)
93 c28(i)=z28(i)+z35(i)
94 20 CONTINUE
95 DO 30 i=lft,llt
96 sx(i)=a17(i)+a28(i)
97 sy(i)=b17(i)+b28(i)
98 sz(i)=c17(i)+c28(i)
99 tx(i)=a17(i)-a28(i)
100 ty(i)=b17(i)-b28(i)
101 tz(i)=c17(i)-c28(i)
102 30 CONTINUE
103
104 DO 100 i=lft,llt
105
106 e3x(i) = ty(i) * rz(i) - tz(i) * ry(i)
107 e3y(i) = tz(i) * rx(i) - tx(i) * rz(i)
108 e3z(i) = tx(i) * ry(i) - ty(i) * rx(i)
109
110 det = sqrt(e3x(i)*e3x(i) + e3y(i)*e3y(i) + e3z(i)*e3z(i))
111 IF ( det/=zero) det = one / det
112 e3x(i) = e3x(i) * det
113 e3y(i) = e3y(i) * det
114 e3z(i) = e3z(i) * det
115
116 c1=sqrt(tx(i)*tx(i)+ty(i)*ty(i)+tz(i)*tz(i))
117 c2=sqrt(rx(i)*rx(i)+ry(i)*ry(i)+rz(i)*rz(i))
118 e1x(i)=tx(i)*c2 +(ry(i) * e3z(i) - rz(i) * e3y(i))*c1
119 e1y(i)=ty(i)*c2 +(rz(i) * e3x(i) - rx(i) * e3z(i))*c1
120 e1z(i)=tz(i)*c2 +(rx(i) * e3y(i) - ry(i) * e3x(i))*c1
121 det = sqrt(e1x(i)*e1x(i) + e1y(i)*e1y(i) + e1z(i)*e1z(i))
122 IF ( det/=zero) det = one / det
123 e1x(i) = e1x(i)*det
124 e1y(i) = e1y(i)*det
125 e1z(i) = e1z(i)*det
126
127 e2x(i) = e3y(i) * e1z(i) - e3z(i) * e1y(i)
128 e2y(i) = e3z(i) * e1x(i) - e3x(i) * e1z(i)
129 e2z(i) = e3x(i) * e1y(i) - e3y(i) * e1x(i)
130 100 CONTINUE
131
132 RETURN