# 實作轉換矩陣

November 27, 2021

## 實作 Matrix3D

NumPy 陣列本身可以作為陣列來使用，而且實作了 `@`（Python 3.5 新增，對應的特殊方法是 `__matmul__`），可用於矩陣相乘：

``````>>> import numpy
>>> a = numpy.array([
...     [1, 2],
...     [3, 4]
... ])
>>> b = numpy.array([
...     [5, 6],
...     [7, 8]
... ])
>>> a @ b # 矩陣相乘
array([[19, 22],
[43, 50]])
>>> a * b # 逐元素相乘
array([[ 5, 12],
[21, 32]])
>>>
``````

``````import numpy

from math import cos, sin, radians

class Matrix3D:
def __init__(self, m):
# self.wrapped 是包裹的 numpy.ndarray 實例
if isinstance(m, numpy.ndarray):
self.wrapped = m
else:
self.wrapped = numpy.array(m)

# Post-Multiplication (Right-Multiplication)
def __matmul__(self, that):
return Matrix3D(self.wrapped @ that.wrapped)

def transform(self, v):
vt = (v.x, v.y, v.z, 1) if isinstance(v, Vector) else v + (1,)
return tuple((self.wrapped @ vt)[:-1])
``````

``````[
row1
row2
row3
]
``````

``````[
[1, 0, tx],
[0, 1, ty],
[0, 0, 1]
]
``````
``````# 單位矩陣
_identity = [
[1, 0, 0, 0],
[0, 1, 0, 0],
[0, 0, 1, 0],
[0, 0, 0, 1]
]

def identity():
return numpy.array(_identity)

# 平移矩陣
def translation(v):
return Matrix3D(_translation(v))

# 繞 x 軸旋轉矩陣
def rotationX(angle):
return Matrix3D(_rotationX(angle))

# 繞 y 軸旋轉矩陣
def rotationY(angle):
return Matrix3D(_rotationY(angle))

# 繞 z 軸旋轉矩陣
def rotationZ(angle):
return Matrix3D(_rotationZ(angle))

# 繞指定軸旋轉矩陣
def rotation(direction, angle):
return Matrix3D(_rotation(direction, angle))

def _translation(v):
vt = (v.x, v.y, v.z) if isinstance(v, Vector) else v
return numpy.array([
[1, 0, 0, vt[0]],
[0, 1, 0, vt[1]],
[0, 0, 1, vt[2]],
[0, 0, 0, 1]
])

def _rotationX(angle):
return numpy.array([
[1, 0, 0, 0],
[0, c, -s, 0],
[0, s, c, 0],
[0, 0, 0, 1]
])

def _rotationY(angle):
return numpy.array([
[c, 0, s, 0],
[0, 1, 0, 0],
[-s, 0, c, 0],
[0, 0, 0, 1]
])

def _rotationZ(angle):
return numpy.array([
[c, -s, 0, 0],
[s, c, 0, 0],
[0, 0, 1, 0],
[0, 0, 0, 1]
])

# 四元數旋轉公式實作
def _rotation(direction, angle):
dir = direction if isinstance(direction, Vector) else Vector(*direction)
axis = dir.normalized()
s = sin(half_a)
x = s * axis.x
y = s * axis.y
z = s * axis.z
w = cos(half_a)
x2 = x + x
y2 = y + y
z2 = z + z
xx = x * x2
yx = y * x2
yy = y * y2
zx = z * x2
zy = z * y2
zz = z * z2
wx = w * x2
wy = w * y2
wz = w * z2
return numpy.array([
[1 - yy - zz, yx - wz, zx + wy, 0],
[yx + wz, 1 - xx - zz, zy - wx, 0],
[zx - wy, zy + wx, 1 - xx - yy, 0],
[0, 0, 0, 1]
])
``````

## 旋轉四面體

``````import numpy
from math import cos, sin, radians
from cadquery import Vector, Edge, Wire, Solid, Shell, Face

class Matrix3D:
def __init__(self, m):
if isinstance(m, numpy.ndarray):
self.wrapped = m
else:
self.wrapped = numpy.array(m)

# Post-Multiplication (Right-Multiplication)
def __matmul__(self, that):
return Matrix3D(self.wrapped @ that.wrapped)

def transform(self, v):
vt = (v.x, v.y, v.z, 1) if isinstance(v, Vector) else v + (1,)
return tuple((self.wrapped @ vt)[:-1])

# 單位矩陣
_identity = [
[1, 0, 0, 0],
[0, 1, 0, 0],
[0, 0, 1, 0],
[0, 0, 0, 1]
]

def identity():
return numpy.array(_identity)

# 平移矩陣
def translation(v):
return Matrix3D(_translation(v))

# 繞 x 軸旋轉矩陣
def rotationX(angle):
return Matrix3D(_rotationX(angle))

# 繞 y 軸旋轉矩陣
def rotationY(angle):
return Matrix3D(_rotationY(angle))

# 繞 z 軸旋轉矩陣
def rotationZ(angle):
return Matrix3D(_rotationZ(angle))

# 繞指定軸旋轉矩陣
def rotation(direction, angle):
return Matrix3D(_rotation(direction, angle))

def _translation(v):
vt = (v.x, v.y, v.z) if isinstance(v, Vector) else v
return numpy.array([
[1, 0, 0, vt[0]],
[0, 1, 0, vt[1]],
[0, 0, 1, vt[2]],
[0, 0, 0, 1]
])

def _rotationX(angle):
return numpy.array([
[1, 0, 0, 0],
[0, c, -s, 0],
[0, s, c, 0],
[0, 0, 0, 1]
])

def _rotationY(angle):
return numpy.array([
[c, 0, s, 0],
[0, 1, 0, 0],
[-s, 0, c, 0],
[0, 0, 0, 1]
])

def _rotationZ(angle):
return numpy.array([
[c, -s, 0, 0],
[s, c, 0, 0],
[0, 0, 1, 0],
[0, 0, 0, 1]
])

# 四元數旋轉公式實作
def _rotation(direction, angle):
dir = direction if isinstance(direction, Vector) else Vector(*direction)
axis = dir.normalized()
s = sin(half_a)
x = s * axis.x
y = s * axis.y
z = s * axis.z
w = cos(half_a)
x2 = x + x
y2 = y + y
z2 = z + z
xx = x * x2
yx = y * x2
yy = y * y2
zx = z * x2
zy = z * y2
zz = z * z2
wx = w * x2
wy = w * y2
wz = w * z2
return numpy.array([
[1 - yy - zz, yx - wz, zx + wy, 0],
[yx + wz, 1 - xx - zz, zy - wx, 0],
[zx - wy, zy + wx, 1 - xx - yy, 0],
[0, 0, 0, 1]
])

def polyhedron(points, faces):
def _edges(vectors, face_indices):
leng_vertices = len(face_indices)
return (
Edge.makeLine(
vectors[face_indices[i]],
vectors[face_indices[(i + 1) % leng_vertices]]
)
for i in range(leng_vertices)
)

vectors = [Vector(*p) for p in points]

return Solid.makeSolid(
Shell.makeShell(
Face.makeFromWires(
Wire.assembleEdges(
_edges(vectors, face_indices)
)
)
for face_indices in faces
)
)

points = (
(5, -5, -5), (-5, 5, -5), (5, 5, 5), (-5, -5, 5)
)

faces = (
(0, 1, 2), (0, 3, 1), (1, 3, 2), (0, 2, 3)
)

n = 6
a_step = 360 / n
for i in range(n):
# 建立轉換矩陣
m = rotationZ(i * a_step) @ translation((15, 0, 0)) @ rotationZ(45)
# 將 points 逐一轉換並建立多面體
solid = polyhedron([m.transform(p) for p in points], faces)
show_object(solid)
``````