行列計算

4次行列の積

4次正方行列m1とm2の積mのi行j列目の計算は以下のようになります。

m[i][j] = m1[i][0]*m2[0][j] + m1[i][1]*m2[1][j] + m1[i][2]*m2[2][j] + m1[i][3]*m2[3][j]

上のiとjを0から3までくりかえします。

頂点座標p1と4次正方行列mの積p2のi列目の計算は以下のようになります。
ただし、p1は[x, y, z, 1]の形をしている斉次座標です。

p2[i] = p1[0]*m[0][i] + p1[1]*m[1][i] + p1[2]*m[2][i] + p1[3]*m[3][i]

上のiを0から3までくりかえします。

行列の積を計算する関数

以下の関数を作成してください。

matmul(m1, m2) 　　

引数

m1, m2 : 4次正方行列

返り値 : 4次正方行列

m1とm2の積

pmatmul(p, m) 　　

引数

p : 斉次座標

m : 4次正方行列

返り値 : 斉次座標

pとmの積

(計算例)

p = [5, 6, 7, 1]
m1 = [[1, 2, 3, 4],
      [5, 6, 7, 8],
      [9, 10, 11, 12],
      [13, 14, 15, 16]]
m2 = [[21,22, 23, 24],
      [25, 26, 27, 28],
      [29, 30, 31, 32],
      [33, 34, 35, 36]]
m = matmul(m1, m2)
print m
print pmatmul(p, m)

とするとm1とm2の積、pとの積が以下のように出力されるようにします。

線形変換

移動、回転、スケールを行う行列は以下のように定義します。

移動

(tx, ty, tz)移動する行列

t = [[1, 0, 0, 0],
     [0, 1, 0, 0],
     [0, 0, 1, 0],
     [tx, ty, tz, 1]]

X軸回転

X軸方向にrラジアン回転する行列

import math
s = math.sin(r)
c = math.cos(r)
rx = [[1, 0, 0, 0],
      [0, c, s, 0],
      [0, s, c, 0],
      [0, 0, 0, 1]]

Y軸回転

Y軸方向にrラジアン回転する行列

import math
s = math.sin(r)
c = math.cos(r)
rx = [[c, 0, -s, 0],
      [0, 1, 0, 0],
      [s, 0, c, 0],
      [0, 0, 0, 1]]

Z軸回転

Z軸方向にrラジアン回転する行列

import math
s = math.sin(r)
c = math.cos(r)
rx = [[c, s, 0, 0],
      [-s, c, 0, 0],
      [0, 0, 1, 0],
      [0, 0, 0, 1]]

スケール

(sx, sy, sz)スケールする行列

s = [[sx, 0, 0, 0],
     [0, sy, 0, 0],
     [0, 0, sz, 0],
     [0, 0, 0, 1]]

計算例

頂点p=[1, 2, 3]をX軸まわりに60度(3.14/3ラジアン)回転してから(2, 3, 4)だけスケールをして(5, 6, 7)移動するには以下のように計算します。

import math

p1 = [1, 2, 3, 1]
s = math.sin(3.14/3)
c = math.cos(3.14/3)
rx = [[1, 0, 0, 0],
      [0, c, s, 0],
      [0, s, c, 0],
      [0, 0, 0, 1]]
s = [[2, 0, 0, 0],
     [0, 3, 0, 0],
     [0, 0, 4, 0],
     [0, 0, 0, 1]]
t = [[1, 0, 0, 0],
     [0, 1, 0, 0],
     [0, 0, 1, 0],
     [5, 6, 7, 1]]

# p*rx*s*t
p2 = pmatmul(p1, rx)
p2 = pmatmul(p2, s)
p2 = pmatmul(p2, t)
print p2

# p*(rx*s*t)
m = matmul(rx, s)
m = matmul(m, t)
p3 = pmatmul(p1, m)
print p3

# p2とp3は同じ値が表示される

結果は以下のように出力されます。

以下のクラスをmatrix.pyというファイル名で作成してください。

class Coodinate4(object):

斉次座標を表すクラス

class Matrix4(object):

4次行列を表すクラス
以下のクラスの親クラスになります。

class TranslateM(Matrix4):

移動行列を表すクラス

class RotateXM(Matrix4):

回転Xを行う行列を表すクラス

class RotateYM(Matrix4):

回転Yを行う行列を表すクラス

class RotateZM(Matrix4):

回転Zを行う行列を表すクラス

class ScaleM(Matrix4):

スケールを行う行列を表すクラス

(計算例)

p = Coodinate4(1, 2, 0)
rx = RotateXM(3.14/2)
ry = RotateYM(3.14/2)
rz = RotateZM(3.14/2)
print p*rx # -> [1, 0.0015926534214665267, 1.9999993658636692, 1]
print p*ry # -> [0.00079632671073326335, 2, -0.99999968293183461, 1]
print p*rz # -> [-1.9992030391529361, 1.001592336353301, 0, 1]

p = Coodinate4(1, 2, 3)
rx = RotateXM(3.14/3)
ry = RotateYM(3.14/3)
rz = RotateZM(3.14/3)
s = ScaleM(2, 3, 4)
t = TranslateM(5, 6, 7)
print p*rx*ry*rz*s*t # -> [2.0720547618139817, 19.971655947266981, 10.008702646121966, 1.0]
print p*(rx*ry*rz*s*t) # -> [2.0720547618139817, 19.971655947266981, 10.008702646121966, 1.0]