したかみ ぶろぐ

Unity成分多め

3D数学の復習と実践(Quaternion)

始めに

過去記事

shitakami.hatenablog.com

shitakami.hatenablog.com

shitakami.hatenablog.com

当初の目標であるQuaternionを扱っていきます。恐らく今回でラストの予定です。



虚数複素数

Quaternionでは複素数を扱うので簡単に解説します。

複素数では虚数を扱います。虚数-1の平方根です。


i^2 = -1


この虚数を利用して、複素数を表現します。2つの実数 a, b を使って次のように表現します。a を z の実数部、b を虚部と言います。


z = a + bi


複素数複素平面で表現が可能で横軸が実軸、縦軸が虚軸です。

f:id:vxd-naoshi-19961205-maro:20191204172529p:plain


原点から複素数を結んだ線と実軸からなる角度を θ とすると次のように書き換えることが出来ます。 ( r は絶対値、大きさ、ノルムと言います。今回は大きさで話していきます)


r = \sqrt{a^2 + b^2} \\
z = r(\cos \theta + i \sin \theta)


複素数はとても便利なことに、絶対値が1の点を原点を軸に θ' だけ回転したいとすると次のように書くことが出来ます。

f:id:vxd-naoshi-19961205-maro:20191204173300p:plain


(\cos \theta + i\sin \theta) \cdot (\cos \theta' + i\sin \theta') \\
= (\cos \theta \cdot \cos \theta'a - \sin \theta \cdot \sin \theta') + 
i(\sin \theta \cdot \cos \theta' + \cos \theta \cdot \sin \theta') \\
= \cos (\theta + \theta') + i \sin (\theta + \theta')



Quaternion

では、次に3次元での回転を扱っていきます。因みにですが、Quaternionは四元数と訳されます。

四元数では3つの虚数 i, j, k が登場します。


i^2 = j^2 = k^2 = -1 \\
ij = k  \quad   ji = -k \\
jk = i  \quad  kj = -i \\
ki = j  \quad  ik = -j


この3つの虚数を使用してQuaternionは次のように表現されます。


\mathbf{q} =
\begin{bmatrix}
w & \mathbf{n}
\end{bmatrix} 
= 
\begin{bmatrix}
w & 
\begin{pmatrix}
x & y & z
\end{pmatrix}
\end{bmatrix} \\
= w + xi + yj + zk


ここからQuaternionの計算について簡単にまとめていきます。


符号反転


-\mathbf{q} =
\begin{bmatrix}
-w & -\mathbf{n}
\end{bmatrix} \\\
=
\begin{bmatrix}
-w &
\begin{pmatrix}
-x & -y & -z
\end{pmatrix}
\end{bmatrix}


大きさ


||\mathbf{q}|| = \sqrt{w^2 + x^2 + y^2 + z^2} \\\
= \sqrt{w^2 + ||n||}


共役


\mathbf{q^*} =
\begin{bmatrix}
w &
\mathbf{-n}
\end{bmatrix} \\\
=
\begin{bmatrix}
w &
\begin{pmatrix}
-x & -y & -z
\end{pmatrix}
\end{bmatrix}


逆数


\mathbf{q^{-1}} =  \displaystyle\frac{\mathbf{q^*}}{||\mathbf{q}||}


乗算(外積


\mathbf{q_{1}} = \begin{bmatrix}
w_{1} & \mathbf{n_{1}}
\end{bmatrix} 
=
\begin{bmatrix}
w_{1} & 
\begin{pmatrix}
x_{1} & y_{1} & z_{1}
\end{pmatrix}
\end{bmatrix}
\\\
\mathbf{q_{2}} = \begin{bmatrix}
w_{2} & \mathbf{n_{2}}
\end{bmatrix} 
=
\begin{bmatrix}
w_{2} & 
\begin{pmatrix}
x_{2} & y_{2} & z_{2}
\end{pmatrix}
\end{bmatrix}

とすると


\mathbf{q_{1}}\mathbf{q_{2}}
=
(w_{1} + x_{1}i + y_{1}j + z_{1}k)(w_{2} + x_{2} + y_{2}j + z_{2}k) \\\
=w_{1}w_{2} - x_{1}x_{2} - y_{1}y_{2} - z_{1}z_{2} \\\
+ (w_{1}x_{2} + w_{2}x_{1} + y_{1}z_{2} - y_{2}z_{1})i  \\\
+ (w_{1}y_{2} + w_{2}y_{1} + z_{1}x_{2} - z_{2}x_{1})j  \\\
+ (w_{1}z_{2} + w_{2}z_{1} + x_{1}y_{2} - x_{2}y_{1})k \\\
= 
\begin{bmatrix}
w_{1}w_{2} - x_{1}x_{2} - y_{1}y_{2} - z_{1}z_{2} \\
\begin{pmatrix}
(w_{1}x_{2} + w_{2}x_{1} + y_{1}z_{2} - y_{2}z_{1}) \\
(w_{1}y_{2} + w_{2}y_{1} + z_{1}x_{2} - z_{2}x_{1}) \\
(w_{1}z_{2} + w_{2}z_{1} + x_{1}y_{2} - x_{2}y_{1}) \\
\end{pmatrix}
\end{bmatrix}


また、次のように表現することも出来ます。


\mathbf{q_{1}}\mathbf{q_{2}}
= \begin{bmatrix}
w_{1} & \mathbf{n_{1}}
\end{bmatrix} 
\begin{bmatrix}
w_{2} & \mathbf{n_{2}}
\end{bmatrix} \\\
=
\begin{bmatrix}
w_{1}w_{2} - \mathbf{n_{1}} \cdot \mathbf{n_{2}} &
w_{1}\mathbf{n_{2}} + w_{2}\mathbf{n_{1}} + \mathbf{n_{1}} \times \mathbf{n_{2}}
\end{bmatrix}


Quaternionの乗算は結合法則は成り立ちますが、交換法則は成り立ちません。


(\mathbf{a}\mathbf{b})\mathbf{c} = \mathbf{a}(\mathbf{b}\mathbf{c}) \\\
\mathbf{a}\mathbf{b} \neq \mathbf{b}\mathbf{a}


2つのQuaternionの積の大きさは次のようになります。 これは2次元の複素平面と同様です。


||\mathbf{q_{1}}\mathbf{q_{2}}|| = ||\mathbf{q_{1}}|| \hspace{3pt} ||\mathbf{q_{2}}||


2つのQuaternionの積の逆数は順序を逆にした逆数の積と等しくなります。


(\mathbf{a}\mathbf{b})^{-1} = \mathbf{b}^{-1}\mathbf{a}^{-1}


内積

Quaternionの内積はベクトルと同様にスカラーとなります。


\mathbf{a} \cdot \mathbf{b} = 
\begin{bmatrix} w_{1} & \begin{pmatrix} x_{1} & y_{1} & z_{1} \end{pmatrix} \end{bmatrix} 
\cdot 
\begin{bmatrix} w_{2} & \begin{pmatrix} x_{2} & y_{2} & z_{2} \end{pmatrix} \end{bmatrix} \\\
=
w_{1}w_{2} + x_{1}x_{2} + y_{1}y_{2} + z_{1}z_{2}


また、2つのQuaternionの大きさがそれぞれ1の場合、値域が-1から1になります。(たまに使える)

 
||\mathbf{a}|| = ||\mathbf{b}|| = 1 \\\
-1 \leqq ||\mathbf{a} \cdot \mathbf{b}||  \leqq 1



Quaternionを使用した回転表現

ここからQuaternionで回転を表現してきます。軸を  \mathbf{n} 、回転させる量を  \theta とすると次のようになります。 また、軸の大きさは1とします。

 
\mathbf{q} =
\begin{bmatrix}
\cos \frac{\theta}{2} & \sin \frac{\theta}{2} \mathbf{n}
\end{bmatrix} \\\
= \begin{bmatrix}
\cos \frac{\theta}{2} &
\begin{pmatrix}
\sin \frac{\theta}{2}n_{x} &
\sin \frac{\theta}{2}n_{y} &
\sin \frac{\theta}{2}n_{z}
\end{pmatrix}
\end{bmatrix} 
\\\
\\\
||\mathbf{n}|| = 1


3次元上の点P(x, y, z)をQuaternion  \mathbf{p} に置き換えて、軸  \mathbf{n} に対して  \theta だけ回転させて得られる点 p' は次のようになります。

 
\mathbf{p} = 
\begin{bmatrix}
0 & \mathbf{v}
\end{bmatrix}
\\\
= 
\begin{bmatrix}
0 & 
\begin{pmatrix}
x & y & z
\end{pmatrix}
\end{bmatrix} 
=
xi + yj + zk
\\\
\\\
\mathbf{p'} =
\mathbf{q}\mathbf{p}\mathbf{q^{-1}}


では、計算して回転を表現しているか確かめていきます。

 
\mathbf{q}\mathbf{p} 
= 
\begin{bmatrix}
\cos \frac{\theta}{2} & \sin \frac{\theta}{2} \mathbf{n}
\end{bmatrix} 
\begin{bmatrix}
0 & \mathbf{v}
\end{bmatrix}
\\\
=
\begin{bmatrix}
0 - \sin \frac {\theta}{2}(\mathbf{n} \cdot \mathbf{v}) & 
\cos \frac {\theta}{2} \mathbf{v} + 0 + \sin \frac{\theta}{2} (\mathbf{n} \times \mathbf{v})
\end{bmatrix} \\\
=
\begin{bmatrix}
 - \sin \frac {\theta}{2}(\mathbf{n} \cdot \mathbf{v}) & 
\cos \frac {\theta}{2} \mathbf{v} + \sin \frac{\theta}{2} (\mathbf{n} \times \mathbf{v})
\end{bmatrix}


ここからの計算は量が増えるので実部と虚部に分けて計算していきます。

 
\mathbf{p'} =
\begin{bmatrix} w' & \mathbf{v'}
\end{bmatrix}
\\\
= 
\begin{bmatrix}
 - \sin \frac {\theta}{2}(\mathbf{n} \cdot \mathbf{v}) & 
\cos \frac {\theta}{2} \mathbf{v} + \sin \frac{\theta}{2} (\mathbf{n} \times \mathbf{v})
\end{bmatrix}
\begin{bmatrix}
\cos \frac{\theta}{2} & -\sin \frac{\theta}{2} \mathbf{n}
\end{bmatrix}


 w' の計算

 
w' = 
- \sin \frac{\theta}{2} \cos \frac{\theta}{2} (\mathbf{n} \cdot \mathbf{v}) 
-  (\cos \frac{\theta}{2} \mathbf{v} + \sin \frac{\theta}{2} (\mathbf{n} \times \mathbf{v})) 
\cdot (-\sin \frac{\theta}{2} \mathbf{n})
\\\
=
-\sin \frac{\theta}{2} \cos \frac{\theta}{2} (\mathbf{n} \cdot \mathbf{v}) +
\sin \frac{\theta}{2} \cos \frac{\theta}{2} (\mathbf{n} \cdot \mathbf{v}) +
\sin ^2 \frac{\theta}{2} (\mathbf{n} \times \mathbf{v}) \cdot \mathbf{n}
\\\
= \sin ^2 \frac{\theta}{2} (\mathbf{n} \times \mathbf{v}) \cdot \mathbf{n}


ここで  \mathbf{n} \times \mathbf{v} \mathbf{n} は直交するため

 
 (\mathbf{n} \times \mathbf{v}) \cdot \mathbf{n} = 0


よって、実部の値が 0 になります。

 
w' =  \sin ^2 \frac{\theta}{2}  (\mathbf{n} \times \mathbf{v}) \cdot \mathbf{n} = 0


 \mathbf{v'}の計算

 
\mathbf{v'} =
\sin^2 \frac{\theta}{2} (\mathbf{n} \cdot \mathbf{v}) \mathbf{n} 
+ \cos \frac{\theta}{2}(\cos \frac{\theta}{2} \mathbf{v} + \sin \frac{\theta}{2}(\mathbf{n} \times \mathbf{v}))
\\\
+ (\cos \frac{\theta}{2}\mathbf{v} + \sin \frac{\theta}{2}(\mathbf{n} \times 
\mathbf{v})) \times (-\sin \frac{\theta}{2}\mathbf{n})
\\\
\\\
= 
\sin^2 \frac{\theta}{2} (\mathbf{n} \cdot \mathbf{v}) \mathbf{n}
+ \cos^2 \frac{\theta}{2} \mathbf{v} +\sin \frac{\theta}{2} \cos \frac{\theta}{2}(\mathbf{n} \times \mathbf{v}) 
\\\
- \sin \frac{\theta}{2} \cos \frac{\theta}{2}(\mathbf{v} \times \mathbf{n}) 
- \sin^2 \frac{\theta}{2}(\mathbf{n} \times \mathbf{v}) \times \mathbf{n}


ここで、外積の順序を入れ替え -1 を掛けたものが元の外積と等しいことから次のように計算できる。

 
\mathbf{v} \times \mathbf{n} = - (\mathbf{n} \times \mathbf{v})
\\\
\\\
\\\
\mathbf{v'} =
\sin^2 \frac{\theta}{2} (\mathbf{n} \cdot \mathbf{v}) \mathbf{n}
+ \cos^2 \frac{\theta}{2} \mathbf{v} +\sin \frac{\theta}{2} \cos \frac{\theta}{2}(\mathbf{n} \times \mathbf{v}) 
\\\
+ \sin \frac{\theta}{2} \cos \frac{\theta}{2}(\mathbf{n} \times \mathbf{v}) 
- \sin^2 \frac{\theta}{2}(\mathbf{n} \times \mathbf{v}) \times \mathbf{n} 
\\\
\\\
= \sin^2 \frac{\theta}{2} (\mathbf{n} \cdot \mathbf{v}) \mathbf{n}
+ \cos^2 \frac{\theta}{2} \mathbf{v} +2\sin \frac{\theta}{2} \cos \frac{\theta}{2}(\mathbf{n} \times \mathbf{v}) 
- \sin^2 \frac{\theta}{2}(\mathbf{n} \times \mathbf{v}) \times \mathbf{n}


ここでベクトルの三重積を使用します。wikipedia等にも解説されているので知らない方は参考にして下さい。

 
(\mathbf{n} \times \mathbf{v}) \times \mathbf{n} 
= (\mathbf{n} \cdot \mathbf{n}) \mathbf{v} - (\mathbf{v} \cdot \mathbf{n}) \mathbf{n}
\\\
= ||\mathbf{n}|| \mathbf{v} - (\mathbf{v} \cdot \mathbf{n}) \mathbf{n}
\\\
= \mathbf{v} - (\mathbf{v} \cdot \mathbf{n}) \mathbf{n}

\\\
\\\
\\\
\mathbf{v'} 
= \sin^2 \frac{\theta}{2} (\mathbf{n} \cdot \mathbf{v}) \mathbf{n}
+ \cos^2 \frac{\theta}{2} \mathbf{v} +2\sin \frac{\theta}{2} \cos \frac{\theta}{2}(\mathbf{n} \times \mathbf{v}) 
- \sin^2 \frac{\theta}{2}( \mathbf{v} - (\mathbf{v} \cdot \mathbf{n}) \mathbf{n})
\\\
= \sin^2 \frac{\theta}{2} (\mathbf{n} \cdot \mathbf{v}) \mathbf{n}
+ \cos^2 \frac{\theta}{2} \mathbf{v} +2\sin \frac{\theta}{2} \cos \frac{\theta}{2}(\mathbf{n} \times \mathbf{v}) 
- \sin^2 \frac{\theta}{2}\mathbf{v} + \sin^2 \frac{\theta}{2}(\mathbf{v} \cdot \mathbf{n}) \mathbf{n}
\\\
=
2\sin^2 \frac{\theta}{2} (\mathbf{n} \cdot \mathbf{v}) \mathbf{n}
+ (\cos^2 \frac{\theta}{2} - \sin^2 \frac{\theta}{2}) \mathbf{v}
+ 2\sin \frac{\theta}{2} \cos \frac{\theta}{2}(\mathbf{n} \times \mathbf{v})


ここで三角関数の二倍角の公式より  \mathbf{v'} は次のようにまとまります。

 
2\sin^2 \frac{\theta}{2} = 1 - \cos \theta \\\
\cos^2 \frac{\theta}{2} - \sin^2 \frac{\theta}{2} = \cos \theta \\\
2\sin \frac{\theta}{2} \cos \frac{\theta}{2} = \sin \theta

\\\
\\\
\\\
\mathbf{v'} 
= (1 - \cos \theta) (\mathbf{n} \cdot \mathbf{v}) \mathbf{n}
+ \cos \theta \mathbf{v}
+ \sin \theta (\mathbf{n} \times \mathbf{v}) 
\\\
= \cos \theta \mathbf{v} + (1 - \cos \theta) (\mathbf{n} \cdot \mathbf{v}) \mathbf{n} + \sin \theta (\mathbf{n} \times \mathbf{v})



以上から求められたQuaternionは次のようになります。

 
\mathbf{p'} =
\begin{bmatrix} 0 & \cos \theta \mathbf{v} + (1 - \cos \theta) (\mathbf{n} \cdot \mathbf{v}) \mathbf{n} + \sin \theta (\mathbf{n} \times \mathbf{v}) 
\end{bmatrix}


ここまで整理出来ました。そして虚部の式は前回解説したロドリゲスの回転行列と一致します。


\vec{v'} 
= \cos \theta \vec{v} + (1 - \cos \theta)(\vec{v} \cdot \vec{n})\vec{n} + \sin \theta (\vec{n} \times \vec{r})


つまり、 \mathbf{p'} =\mathbf{q}\mathbf{p}\mathbf{q^{-1}} は回転を表しています。

このQuaternionの計算で求まった  \mathbf{v'} は座標  \mathbf{v} を軸  \mathbf{n}  \theta 回転した座標となります。



Quaternionの利点・欠点

Quaternionはオイラー角に比べて以下の利点があります。

しかし、欠点として次のものがあります。

  • 私たちにとって理解しずらい。
  • オイラー角と比べると変数が一つ増える



最後に

やっとQuaternionについてまとめることが出来ました。Quaternionが回転を表していることをどうしても理解したかったのでいい機会となりました。

UnityでもQuaternionクラスが実装されており、とても便利なのでオイラー角で求めるのが難しい場合は使うべきものかもしれません。

ちなみに、昔にオイラー角での回転とQuaternionでの回転の計算速度を簡単に計測した結果、Quaternionのほうがわずかに高速だとわかりました。


また、ドローンの動きを作成する際にもQuaternionを使用してドローンが傾くように実装しました。



参考

3D数学の復習と実践の始めのブログで載せた参考書に加えて次の参考書を参考に式を導出しました。

実例で学ぶゲーム3D数学

実例で学ぶゲーム3D数学

Unityでわかる! ゲーム数学

Unityでわかる! ゲーム数学

  • 作者:加藤 潔
  • 出版社/メーカー: 翔泳社
  • 発売日: 2018/06/18
  • メディア: 単行本(ソフトカバー)