Quaternionのお勉強

目次 (2023/8/12 ver.0.04)

1番目のタブが選択された!

1.What's New

目次

1.1 Matrix4について

4D matrix とは
  • 三次元を表現するために4*4の行列で表すものが「4D matrix」です。
  • なので、Matrix4クラスのコンストラクタは16個の引数を必要とします。
  • $\Large \displaystyle Matrix4($

    $\Large \displaystyle \ 1, 0, 0, 0$

    $\Large \displaystyle \ 0, 1, 0, 0$

    $\Large \displaystyle \ 0, 0, 1, 0$

    $\Large \displaystyle \ 0, 0, 0, 1$

    $\Large \displaystyle )$

  • 上記のコードのように、右下がりの対角線上にある成分がすべて1で、残りの成分がすべて0であるものを identity matrix (単位行列)と呼びます。
  • identity matrix をtransformに与えても何も変更されません。

4D matrix を使ってScaleしてみる
  • こちらのデモはボタンを押すとx方向に2倍になります。
  • import 'package:flutter/material.dart';
    
    void main() {
      runApp(
        MaterialApp(
          debugShowCheckedModeBanner: false,
          home: MainPage(),
        ),
      );
    }
    
    class MainPage extends StatefulWidget {
      @override
      _MainPageState createState() => _MainPageState();
    }
    
    class _MainPageState extends State {
    
      double x = 1;
    
      @override
      Widget build(BuildContext context) {
        return Scaffold(
          floatingActionButton: FloatingActionButton(
            onPressed: () {
              setState(() {
                x = 2;
              });
            },
            child: Icon(Icons.add),
          ),
          body: Center(
            child: Transform(
              transform: Matrix4(
                x, 0, 0, 0,
                0, 1, 0, 0,
                0, 0, 1, 0,
                0, 0, 0, 1
              ),
              child: Container(
                width: 100,
                height: 100,
                color: Colors.indigo[500],
              ),
            ),
          ),
        );
      }
    }
    
  • このように行列の値を変化させることでScaleしたりtranslateしたりできます。
  • しかし行列を直接いじるためには知識が必要になってしまうので(Rotationはsin, cosもでてくる)、簡単に扱えるようにMatrix4クラスにはファクトリが定義されています。
いろんなファクトリ
大きさを変える
    Matrix4.diagonal3Values(double x, double y, double z)
    
回転させる
    Matrix4.rotationX(double radians)
    Matrix4.rotationY(double radians)
    Matrix4.rotationZ(double radians)
    
  • X, Y, Zは回転軸を表しています。
  • 引数に渡すのは回転させたい角度ではなくてラジアンで指定します。
  • radian = angle * (pi / 180)
    
  • たとえば90度回転させたい場合は、
  • radian = 90 * (pi / 180) = 0.5 * pi
    
  • となります。
ななめにする ゆがめる
  • 説明が難しいのでデモを見てください。

1.2 使用する記号等

機体固定座標系

2. Quaternionの概要

引用先:「モーションにおける3次元回転 #4 〜オイラーの公式からロドリゲスの式へ〜

目次

2.1 Quaternionについて思うこと

1.説明
本記事ではQuaternionについて総整理します。Quaternion (四元数) は、
  • 航空機
  • 宇宙機
  • ロボティクス
  • 3D ゲーム
  • 3D コンピュータグラフィックス
などなど、幅広い分野で活用されているトピックです。
航空宇宙業界と 3D ゲーム業界との両方で盛んに用いられているという点でとても興味深いです。
特に最近では Unity の流行により、職業ゲームプログラマだけでなく、趣味で 3D ゲーム制作を楽しむ方にとっても、Quaternionについて理解したいという方が増えているのを感じます。
本記事では、Quaternionについて必要な事柄を総整理するとともに、Quaternionをワカッタ感覚で納得することを目標にします。
2.フライトシュミレータ
Unityを使って鳥人間コンテストの滑空機部門を想定したフライトシミュレーターを作成する場合、 以下のような計算が必要になるらしい。
  • 質量・重心位置・慣性モーメント
  • 迎角と横滑り角 (α, β)
  • 突風成分
  • 地面効果
  • 縦のつり合い
  • 縦の空気力 X,Z
  • 縦のモーメント M
  • 横・方向のつり合い
  • 横・方向の空気力 Y
  • 横・方向のモーメント L,N

2.2 オイラーの公式からロドリゲスの式へ

オイラーの定理
オイラーの定理は、3次元において物体の回転を
「ひとつの回転軸と回転角度で表すことができる」
ということを示している。

図1: オイラーの定理

たとえば,任意のベクトル$a$を回転行列$R$によって回転させたベクトル$a'$は,
$Ra=a'$
と書くことができ、図1の点Pは回転 $R$ によって点P'に移動する.

これに対して,唯一回転軸 $n$ に対しては
$Rn=n$
が成り立ち,図1の点Qは同じ位置にとどまる.
これは回転 $R$ によって「不変なベクトル」 $n$ が存在し,それが回転軸であることに相当する.

なお, $n$ を単位ベクトルとし,
$Rn= \lambda n$
で記述した直交行列 $R$ に対する固有値$\lambda$は, 1, $e^{i \psi}$, $e^{-i \psi}$ である。
このうちの1である固有値は,回転を行っても不変を意味することから,唯一の回転軸の存在理由となっている.
ロドリゲスの式
物体の回転は,オイラーの定理から「ひとつの回転軸 $n$ と回転角度 $\psi$ で表せる」ので,回転行列は
$R = E + [n \times ] \sin \psi + [n \times ]^2(1- \cos \psi) $  (1)
のように表せる.これがロドリゲスの式(Rodrigues formula)として知られている.
ここで, $E$ は $3\times3$ の単位行列, $[n \times]$ は外積を表す歪対称行列である(次節参照).この式の導出については最後に示す.
ロドリゲスの式は幾何学条件を使用し,他にも多様な表現が可能だが,ここではこの式の幾何学的意味を考える.
その前に,外積と歪対称行列について次で説明を行う.
外積
ベクトル $a$ と $b$ の外積(outer product)はベクトルであることから $a \times b$ ベクトル積(vector product)とも呼ばれる。
一般に記号で $a \times b$ のように表現されるためクロス積(cross product)とも呼ばれる.
ベクトル $a$ と $b$ の外積(outer product) $a \times b$ は,図2のように,2つのベクトル $a$, $b$ に対して垂直な方向のベクトルを定める.

図2: 外積の幾何学的意味

具体的な外積 $a \times b$ の幾何学的性質は
1. $a$, $b$ の両方に直交し,始点を一致させて $a$ を $b$ に近づけるように,右ネジを回す際の進行方向を向く.
2. その大きさは, $a$, $b$ の始点を一致させて作る平行四辺形の面積($|a||b|\sin\psi$ )に等しい.

また,2つのベクトル
$a=\begin{bmatrix}a_x\\a_y\\a_z\end{bmatrix}$, $b=\begin{bmatrix}b_x\\b_y\\b_z\end{bmatrix}$
に対して,外積の成分を
$a \times b=\begin{bmatrix}a_y b_z - a_z b_y\\a_z b_x - a_x b_z\\a_x b_y - a_y b_x\end{bmatrix}$
のように書くことができる.
歪対称行列を使用した外積表現
$n=[n_x~n_y~n_z]^T$ に対して,先に述べた $[n \times]$ は以下のように定義される行列
$[n \times] \equiv \begin{bmatrix}0&-n_z&n_y\\n_z&0&-n_x\\-n_y&n_x&0\end{bmatrix}$
である.このような $A^T = -A$ の関係を満たす行列を歪対称行列,または交代行列と呼ぶ.
これによって,ベクトルの外積を
$c=a\times b=[a\times]b$
$\begin{bmatrix}c_x\\c_y\\c_z\end{bmatrix}=\begin{bmatrix}a_y b_z - a_z b_y\\a_z b_x - a_x b_z\\a_x b_y - a_y b_x\end{bmatrix}=\begin{bmatrix}0&-a_z&a_y\\a_z&0&-a_x\\-a_y&a_x&0\end{bmatrix}\begin{bmatrix}b_x\\b_y\\b_z\end{bmatrix}$
のように書くことができる.
このような外積計算の行列による表現は,外積計算を含んだ式も線形表現できるので,今後多用することになる.
ロドリゲスの式の行列表現
ここで,ひとつ注意点を述べると,回転行列はオイラー角でも,単位Quaternion(オイラーパラメータ)などで表現されうことが多いが,ロドリゲスの式でも回転行列を表現できることである.
すなわち回転軸を意味する単位ベクトル $n$ の3変数と回転角度 $\psi$ でも,回転行列を表現できることを示している.
これはQuaternionの性質と似ているのだがその意味は,次の章で説明する.

式(1)を展開すると,
$R_{n}(\psi) =\\ \hspace{5mm}\left[\begin{array}{cc}\cos\psi + n_x^2(1-\cos\psi) & n_x n_y(1-\cos\psi)-n_z \sin \psi \\ n_y n_x(1-\cos\psi)+ n_z \sin \psi & \cos\psi + n_y^2(1-\cos\psi) \\ n_z n_x(1-\cos\psi) - n_y \sin \psi & n_z n_y(1-\cos\psi) + n_x \sin \psi \end{array}\right.\\ ~\\ \hspace{6cm} \left.\begin{array} {c}n_x n_z(1-\cos\psi) + n_y \sin \psi\\ n_y n_z(1-\cos\psi) - n_x \sin \psi\\ \cos\psi + n_z^2(1-\cos\psi)\end{array}\right]$
となる.ここで,回転軸の単位ベクトル $n$ の成分を
$n=\begin{bmatrix}n_x\\n_y\\n_z\end{bmatrix}$
とした.
ロドリゲスの式の逆変換
逆に,回転行列が既知のときに回転軸 $n$ と,回転角度 $\psi$ も計算可能である.
解析では回転行列を定めてから回転軸と回転角度を知りたいというこちらのニーズのほうが高いかもしれない.
前述の直交行列の固有値が1の場合を考えると,それに対応する固有ベクトルは回転軸に相当することを述べた.その時の固有ベクトル, すなわち回転軸を $w=[w_1~w_2~w_3]^T$とすると,
$Rw = w$
が成立する.
また, $R$ は直交行列なので, $R^T = R^{-1}R$ が成立し,両辺に $R^T$ を掛けて
$R^T R w = R^Tw \\ w = R^Tw$
となるので,
$(R - R^T)w = 0$
を得る.また,回転行列 $R$ が以下のように既知のとき,
$R=\begin{bmatrix}r_{11}&r_{12}&r_{13}\\ r_{21}&r_{22}&r_{23}\\ r_{31}&r_{32}&r_{33}\end{bmatrix}$
とすると,
$(R - R^T)w = 0 \\ \begin{bmatrix}0 & r_{12}-r_{21} &r_{13} - r_{31}\\ r_{21} - r_{12} & 0 & r_{23} - r_{32} \\ r_{31}-r_{13} & r_{32}-r_{23} & 0\end{bmatrix} \begin{bmatrix}w_1\\w_2\\w_3\end{bmatrix}= \begin{bmatrix}0\\0\\0 \end{bmatrix}$
となるが, $R - R^T$ は歪対称行列となっていることに気がつく.
そこで,これをベクトル $W$ で記述された $[W \times]$ と書くと,これは外積の演算を意味するので,
$W \equiv \begin{bmatrix} r_{32}-r_{23} \\ r_{13} - r_{31} \\ r_{21} - r_{12}\end{bmatrix} \\ R - R^T = \begin{bmatrix}0 & r_{12}-r_{21} &r_{13} - r_{31}\\ r_{21} - r_{12} & 0 & r_{23} - r_{32} \\ r_{31}-r_{13} & r_{32}-r_{23} & 0\end{bmatrix}=[W \times]$
と書け, $W \times w = 0$ を満たす $w$ は,外積の幾何学的意味からも
$w = W$
しかない.そこで,ベクトル $w$ を
$w\equiv\begin{bmatrix}r_{32} - r_{23}\\ r_{13} - r_{31} \\ r_{21}-n_{12}\end{bmatrix}$
と定義すると,大きさが1の回転軸 $n$ は
$n=\frac{w}{||w||}$
のようにベクトル $w$ の単位ベクトルとなり,回転角度 $\psi$ は
$\psi = \mathrm{atan2} (\sin \psi, \cos \psi)$
となる.ただし,ここで
$\cos \psi = \frac{1}{2}(\mathrm{tr}~ R - 1)\\ = \frac{1}{2}(n_{11} + n_{22} + n_{33}-1)$
であり, $n$, $\cos \psi$ が既知なので, $\sin \psi$ はロドリゲスの式 $R = E + [n \times ] \sin \psi + [n \times ]^2(1- \cos \psi)$
(1−cosψ) から算出する.

また,トレース $\rm{tr}~ R$ は,前述の行列 $R$ の対角要素の和
${\mathrm{tr}}~ R = n_{11} + n_{22} + n_{33}$
である.

角度を計算する $\rm{atan2} (y,x)$ はプログラミング言語などの組み込み関数として広く使用されている2引数の逆正接関数( $\tan()$ の逆関数)である.
引数の順番が y,xy,x の順番であることに注意する.
これを使用することで,角度が四象限にわたる符号も得られ,かつ計算精度も高い.

なお,PythonのNumpy場合のatan2()関数としては,np.arctan2(y, x)を使用する.
これは引数が二つで,arctan(y / x)をラジアンで返す.
この角度は,極座標平面において原点から座標(x, y)へのベクトルがx軸の正の方向となす角度(偏角)であり,戻り値は -\pi−π から \piπ (-180度から180度)の間になる.
計算言語によって,引数の順番,定義域に差があるので注意されたい.
読者は角度 $\psi$ の計算は, $\cos \psi = \frac{1}{2}(\rm{tr}~ R - 1)$ の逆関数で計算すれば良いと考えるかもしれないが, 上記の符号と計算精度の理由から,上記のように計算が複雑化しても,原則 $\cos ^{-1}()$ を角度計算では使用してはいけないことを留意しておくべきである.

2.3 ロドリゲスの式の導出

ベクトルの計算準備
ここでは,次節でロドリゲスの式の証明を行うための,射影と反射影というベクトルの幾何学的計算準備を行う.

図3: ベクトルの射影と直交成分

ベクトル $a$ はベクトル $n$ に平行な成分 $a_{||}$ と,ベクトル $n$ に対して直交する成分 $a_{\perp}$ とに分解することができるため,ベクトル $a$ は
$a = a_{||} + a_{\perp}$   (2)
と書くことができ, $a_{||}$ を $a$ の $n$ への射影(rejection), $a_{\perp}$ を $a$ の $n$ からの反射影(rejection)とよぶ.
射影 $a_{||}$ の長さは,ベクトル $a$, $n$ の内積 $a^T n$ で計算できるので,射影 $a_{||}$ は
$a_{||}= (a^T n)n$
と書け,さらに
$a_{||} = (a^T n)n= (n^T a)n=n(n^T a)\\ = (n n^T) a$
と書き換えることができる.ここで,内積は順序を入れ替えても良いというルールなどに基づいて書き換えている.
一方,反射影 $a_{\perp}$ は式(2)を利用して,
$a_{\perp} = a - a_{||}= a- (n^T n) a = (E - n^T n)a$
と書くことができる.
この射影と反射影を利用して,次にロドリゲスの式を幾何学的に導いていく.
ロドリゲスの式の幾何学的意味
この節では,ロドリゲスの式の幾何学的意味を理解するために,ロドリゲスの式の幾何学的な導出を行う.
内積,外積,ベクトル三重積などの性質を使用して証明を試みるが,この節は読み飛ばしても構わない.

図4: ロドリゲスの式の幾何学的意味

回転行列 $R$ にベクトル $a$ をかけると,ベクトル $a$ を角度 $\psi$ だけ回転して $a'$ となる.
これを式で表すと,
$a'=Ra = a + [n \times ] a \sin \psi + [n \times ]^2 a(1- \cos \psi)$   (2)
となる.そこで,図4におけるベクトル $a'$ と,この式の関係について検証する.
ベクトル $a'$ は点P'を指し示している.点P'は,原点Oから,点P,Hを経由してP'に達するので
$a' = \overrightarrow{OP}' = \overrightarrow{OP}+\overrightarrow{PH}+\overrightarrow{HP'}$
と書くことができる.そこで,これを分解して考える.まず,円の半径は同じ長さであるので,
$||\overrightarrow{QR}||=||\overrightarrow{QP'}||=||\overrightarrow{QP}||=||[n\times] a||=||a_{\perp}||=||(E-n^T n)a||$
が成り立つ.そして,まず,
1: $ \overrightarrow{OP}$ はベクトル $a$
$\overrightarrow{OP}=a$
である.
2: $\overrightarrow{PH}$ は, $\overrightarrow{PH}$ の長さに $\overrightarrow{PH}$ 方向の単位ベクトルをかけたものなので,
 2-1: まず,長さ $||\overrightarrow{PH}||$ を求める.
$||\overrightarrow{PH}||=||\overrightarrow{HP}||=||\overrightarrow{QP}-\overrightarrow{QH}||$ であるため.
 2-1-1: まず, $\overrightarrow{QP}$ を考えると,ベクトル $\overrightarrow{QP}$ は $a$ の $n$ からの反射影 $a_{\perp}a$ であるため, 前節の結果を利用し, $\overrightarrow{QP}=a_{\perp}= (E - n^T n)a$ を得る.
 2-1-2: また,ベクトル \overrightarrow{QH} QH ​ は, \overrightarrow{QP’} QP’ ​ のベクトルの \overrightarrow{QP}=\bm{a}_{\perp} QP ​ =a ⊥ ​ への射影であるので
$\overrightarrow{QH} = \cos \psi ||\overrightarrow{QP'}|| \frac{\overrightarrow{QP}}{||\overrightarrow{QP}||}= \cos \psi ||\overrightarrow{QP}|| \frac{\overrightarrow{QP}}{||\overrightarrow{QP}||}= \cos \psi ~\overrightarrow{QP} \\=\cos\psi ~a_{\perp}=\cos \psi(E-n^Tn) a$
となる.2-1-1と2-1-2の結果を合わせて,
$||\overrightarrow{PH}|| =||\overrightarrow{QP}-\overrightarrow{QH}||\\=||a_{\perp}-\cos \psi a_{\perp}||\\=||(1-\cos \psi)a_{\perp}||$
となる.以上から, $\overrightarrow{PH}$ は,長さ $||\overrightarrow{PH}||$ に $\overrightarrow{PQ}$ 方向の単位ベクトルをかけたものなので,
$\overrightarrow{PH} = -||(1-\cos \psi)a_{\perp}|| ~\overrightarrow{QP} = -||(1-\cos \psi)a_{\perp}|| ~\frac{a_{\perp}}{||a_{\perp}||}\\=-(1 - \cos \psi)a_{\perp}= -(1 - \cos \psi)(E-n^Tn) a$
を得る.
3: 最後に,ベクトル $\overrightarrow{HP'}$ は,長さ $||\overrightarrow{HP'}||=||\overrightarrow{QP'}||\sin \psi = ||~[n\times] a~|| \sin \psi$ に対して, $\overrightarrow{HP'}$ の向きのベクトル $[n\times] a$ の単位ベクトル $\frac{[n\times] a}{||~[n\times] a~||} $
をかけたものであるので,
$\overrightarrow{HP'}=||~[n\times] a~|| \sin \psi \frac{[n\times] a}{||~[n\times] a~||} = \sin \psi ~[n\times] a $
を得る.以上をまとめて,
$a' = \overrightarrow{OP}' = \overrightarrow{OP}+\overrightarrow{PH}+\overrightarrow{HP'}\\=a-(1 - \cos \psi)(E-n^Tn) a+\sin \psi ~[n\times] a$  (3)
を得る.
なお,ベクトルの三重積から
$(c\times b) a = ([c \times ][b \times ]) a = (bc^T )a − (c^T b)a$
が成り立つが,ここで $b$, $c$ それぞれに $n$ を代入すると,
$([n \times ][n \times ])a = (nn^T )a − Ea = (nn^T − E)a$
となることを利用し,式(3)に代入すると,
$a' = \overrightarrow{OP}' = \overrightarrow{OP}+\overrightarrow{PH}+\overrightarrow{HP'}\\=a+(1 - \cos \psi)([n \times ][n \times ]) a+\sin \psi ~[n\times] a\\= a + [n \times ] a \sin \psi + [n \times ]^2 a(1- \cos \psi)$
を得て,これは式(1)のロドリゲスの式と一致する.以上,ロドリゲスの式を幾何学的に導出した.
最後に
ロドリゲスの式の幾何学的意味を説明した.
ここで重要なことは,この式は回転行列 $R$ を回転軸を示す単位ベクトル $n$ と回転角度 $\psi$ だけで表現でき,特にこの式(1)の形では,複素数で2次元の回転運動を記述するオイラーの公式を連想させる形式で書かれていることにある.
式(1)のベクトル $\overrightarrow{OP}=a$ で平面に移動した後は,図4右の円で表される平面内の運動を考えればよく,この式は,なにやら $[n\times]$ が複素数における虚数と似たような役割を果たしそうなことと,3次元の回転を表すQuaternionが2次元の複素数の拡張になっている?ことを予感させてくれる.

3. Quaternionを用いてやりたいこと

目次

3.1 Quaternionを使って何をするか

そもそもQuaternionを使って何がやりたいのでしょうか?大きく分けて以下の 2 つがあるように思います。
  • 3D 空間における「回転」を表す
  • 3D 物体の「姿勢」を表す
2 つと言ったのですが、メインは「回転」の方でしょう。
「姿勢」は単純に「どういう回転でその向きになったか」として捉えることができます。
注意点として、「姿勢」は「方向」よりも多くの情報を含んでいます。
Quaternionを用いた回転を扱うときに、回転させたい対象が「姿勢」なのか「方向」なのかを混乱しないようにすることが肝要です。
方向は「進行方向」や「視線」といったものを表していますが、姿勢はそれに加えて、その方向周りの角度の情報も含みます。
航空機でいえば、機首の向いている方向だけでなく、機首に対して機体がどれだけ傾いているのかも指定してはじめて機体の姿勢が完全に決まります。
数理的には以下のようになっています。
自由度 備考
方向 2 進行方向は三次元ベクトル (3 個のパラメータ) で表せますが、そのベクトルの大きさは意味を持たないので 1 減らして 2 です
姿勢 3 方向 + 1 です
つまり、姿勢は最小 3 個のパラメータで表せるということです。
回転も姿勢と本質的に同じものなので、同様に 3 個のパラメータで表せます。
実際にオイラー角は 3 個のパラメータで姿勢や回転を表します (回転を表す方法はQuaternionだけではなく、オイラー角による方法もメジャーな方法の 1 つです)。
オイラー角は Unity ではインスペクターの Rotation 項目で表示されているやつです。
(しかし transform.rotation は Quaternion 型なので少し紛らわしいですね...この記事でも注意喚起がなされています)

3.2 Quaternionとは

Quaternionは任意軸回転がキーワードとして言及される場面も多いです。
いかなる回転 (三次元空間内) も表現できることが強みです。具体的には以下のように定義します。
三次元空間において
  • 方向ベクトル $(\lambda_x, \lambda_y, \lambda_z)$ を回転軸として
  • 角度 $\theta$ だけ回転する
という「回転」を表すQuaternionは、四次元ベクトル $(\lambda_x \sin⁡{\frac{\theta}{2}}, \lambda_y \sin⁡{\frac{\theta}{2}}, \lambda_z \sin⁡{\frac{\theta}{2}}, \cos{\frac{\theta}{2}})$ で表します。
これらは Unity では、Quaternion 型変数を qua として、
  • qua.x = $\lambda_x \sin⁡{\frac{\theta}{2}}$
  • qua.y = $\lambda_y \sin⁡{\frac{\theta}{2}}$
  • qua.z = $\lambda_z \sin⁡{\frac{\theta}{2}}$
  • qua.w = $\cos{\frac{\theta}{2}}$
という風に対応しています。
(が、Unity の Quaternion を扱うときにこれらの成分を直接触ることはあまりないかもです...スクリプトリファレンスにもQuaternionに熟知していない限りは x, y, z, w の値を直接変更しないでくださいと書いてあります)
また Unity で上記の定義に則ってQuaternionを生成するときは Quaternion.AngleAxis を用います (オイラー角など他の方法でQuaternionを生成する方法は後述します)。


// Y 軸 (上方向) まわりに 30 度回転するのを表すQuaternion
Quaternion.AngleAxis(30, Vector3.up)

ところで上の定義を見ると
  • なんでいきなり三角関数が...
  • $\frac{\theta}{2}$ ってなんだよ... $\theta$ でいじゃん...
と思ってしまいそうです。
その答えは「こう定義すれば具体的な回転計算が超楽だから」なのですが、本記事を読み進めて行くことで少しずつ実感していただけたらと思います。

3.3 Quaternionによる回転操作

Quaternionによる「回転」を定義したからには、実際に「三次元ベクトル」や「姿勢」を回転させる計算方法を導出してあげる必要があります。
三次元ベクトル $\vec{v}$ の回転計算は
  • キャラクターの視線をある方向へと徐々に向けて行く (3D ゲーム)
  • 宇宙機に搭載した望遠鏡を所望の方向へ向ける (航空宇宙)
といった場面で必要になります。
また、姿勢 $p$ の回転計算は、
  • QuaternionとQuaternionの補間 (3D ゲーム)
  • 航空機や宇宙機の三軸姿勢制御 (航空宇宙)
といった場面で必要になります。
結論としては、三次元ベクトル $\vec{v}$ や姿勢Quaternion $p$ に対し、回転 $q$ を実施した結果は以下のように大変シンプルなものになります!!!!!!!
回転させるもの 回転結果 備考
$\vec{v}$ $q \otimes \vec{v} \otimes \bar{q}$ $\bar{q}$ は $q$ の逆回転
姿勢Quaternion $p$ $q \otimes p$
逆に
回転操作をシンプルな計算で実現できるように、Quaternionをあのように定義した
とも言えます。
なお、Quaternion同士のかけ算 $q \otimes p$ や、Quaternionとベクトルとのかけ算 $q \otimes \vec{v}$ はまだ定義していないですが、1-4 で導入します。
また、こうした事柄はQuaternionを自分で実装する場合には必要な知見なのですが、Unity ではいずれも単に
Unityでの回転操作
v2 = q * v   // 三次元ベクトル v の回転
p2 = q * p   // 姿勢 p の回転
と「*」演算子で実現することができます。
$q \otimes \vec{v} \otimes \bar{q}$ な量を q * v と記述するのは少し紛らわしさがありますが、実用上とても簡便です。
演算子の順番は重要で、左右を入れ替えてはダメです。

3.4 なぜQuaternionなのか

三次元空間上での回転や姿勢を表現する方法には代表的なものが 3 つあります。
  • Quaternion
  • オイラー角
  • 回転行列 or DCM
パラメータ数 メモリ消費 計算上の扱いやすさ 備考
Quaternion 4 バランス型で、その他様々なメリットがあります
オイラー角 3 × 具体的な姿勢をイメージしやすく、メモリ的にも優しいですが、計算上の難点を抱えています
回転行列 or DCM 9 ×
どの方法にもメリット・デメリットがあるので、状況に応じて使いこなせばいいのですが、総じてQuaternionはバランスがとれた優秀な方法であることがわかります。

オイラー角はたった 3 個のパラメータで回転や姿勢を表現できるので、メモリ消費は少なくて済むのですが、三角関数を多用することもあって計算時間的には大きくなってしまいます。
それどころかもっと悪いことに、オイラー角の二番目の角度が ±90±90 度になる場合が特異点となって計算上の大きな問題を引き起こすことが知られています (詳しくはジンバルロックを参照)。
端的に言えば、ヨーとロールの区別ができなくなってしまう状態です。
しかしながら特に航空機の世界では
  • そもそもジンバルロックが起こるような状態にならない (ピッチが ±90±90 度になる状態...下図のような墜落時などに発生しそうな相当ヤバイ状態です...戦闘機などでは考慮が必要です)
  • ヨー、ピッチ、ロールという概念がとてもなじみ深い

といった理由によって盛んに用いられています。
また、3D ゲームの世界などにおいても、オイラー角の直感的なイメージしやすさから、UI や API 部分では用いられるケースも多いです。

3.5 Quaternionのかけ算

4.Quaternionの性質

4.1 Quaternionとベクトル

ベクトルの積
  • Quaternionの積を考える前に、ベクトル「$\vec{a}$」と「$\vec{b}$」に以下の新しい積を導入する。
  • $\Large \displaystyle \vec{a}*\vec{b}=\vec{a}\times\vec{b}-\vec{a}\cdot\vec{b}$

  • 「$\vec{a}\times\vec{b}$」は、ベクトルの外積で、「$\vec{a}\cdot\vec{b}$」は、ベクトルの内積である。
  • 以下のように書く場合もあるらしい。
  • $\Large \displaystyle \vec{a}\vec{b}=\vec{a}\times\vec{b}-\vec{a}\cdot\vec{b}$

  • 今、例として「$\vec{a}=\vec{b}=\vec{i}$」とすると、以下のようになる。
  • $\Large \displaystyle \vec{i}*\vec{i}=\vec{i}\times\vec{i}-\vec{i}\cdot\vec{i}=\vec{0}-1=-1$

  • 「$\vec{a}=\vec{i}$」「$\vec{b}=\vec{j}$」とすると、以下のようなる。
  • $\Large \displaystyle \vec{i}*\vec{j}=\vec{i}\times\vec{j}-\vec{i}\cdot\vec{j}=\vec{k}-0=\vec{k}$

  • 「$\vec{j}$」「$\vec{k}$」についても同様の演算を行うと、「$\vec{i}$」「$\vec{j}$」「$\vec{k}$」は、以下のルールに従う。
  • $\Large \displaystyle \vec{i}*\vec{i}=\vec{j}*\vec{j}=\vec{k}*\vec{k}=\vec{i}*\vec{j}*\vec{k}=-1 \ \ \ \ (1)$

    $\Large \displaystyle \vec{i}*\vec{j}=-\vec{j}*\vec{i}=\vec{k} \ \ \ \ (2)$

    $\Large \displaystyle \vec{j}*\vec{k}=-\vec{k}*\vec{j}=\vec{i} \ \ \ \ (2)$

    $\Large \displaystyle \vec{k}*\vec{i}=-\vec{i}*\vec{k}=\vec{j} \ \ \ \ (2)$

  • (1)から、「$\vec{i}$」「$\vec{j}$」「$\vec{k}$」は、虚数単位の「$i^2=-1$」と同じ性質を持つ事が分かる。
  • (2)から、「$\vec{i}$」「$\vec{j}$」「$\vec{k}$」は、ユークリッドの基底ベクトルの演算ルールに一致している事が分かる。
Quaternionの積
  • 次に二つのQuaternionの積を考える。
  • $\Large \displaystyle r=p*q=(p_0 + p_1\vec{i} + p_2\vec{j} + p_3\vec{k})*(q_0 + q_1\vec{i} + q_2\vec{j} + q_3\vec{k})$

    $\Large \displaystyle =p_0 q_0 + p_0 q_1\vec{i} + p_0 q_2\vec{j} + p_0 q_3\vec{k}$

    $\Large \displaystyle +p_1 q_0\vec{i} + p_1 q_1\vec{i}*\vec{i} + p_1 q_2\vec{i}*\vec{j} + p_1 q_3\vec{i}*\vec{k}$

    $\Large \displaystyle +p_2 q_0\vec{j} + p_2 q_1\vec{j}*\vec{i} + p_2 q_2\vec{j}*\vec{j} + p_2 q_3\vec{j}*\vec{k}$

    $\Large \displaystyle +p_3 q_0\vec{k} + p_3 q_1\vec{k}*\vec{i} + p_3 q_2\vec{k}*\vec{j} + p_3 q_3\vec{k}*\vec{k}$

  • 上式に(1)(2)を適用して右辺を整理すると、二つのQuaternionの積「$r$」は、次のようになる。
  • $\Large \displaystyle r=p*q=r_0 + r_1\vec{i} + r_2\vec{j} + r_3\vec{k}$

    $\Large \displaystyle =p_0 q_0 - q_1 q_1 - q_2 q_2 - q_3 q_3$

    $\Large \displaystyle +( p_1 q_0 + p_0 q_1 - p_3 q_2 + p_2 q_3)\vec{i}$

    $\Large \displaystyle +(p_2 q_0 + p_3 q_1 + p_0 q_2 - p_1 q_3)\vec{j}$

    $\Large \displaystyle +(p_3 q_0 - p_2 q_1 + p_1 q_2 + p_0 q_3 )\vec{k}$

  • 両辺の係数を見比べて行列の形にすると、次のようになる。
  • $\Large \displaystyle \begin{bmatrix} r_0 \\ r_1 \\ r_2 \\ r_3 \\ \end{bmatrix} = \begin{bmatrix} p_0 & -p_1 & -p_2 & -p_3\\ p_1 & p_0 & -p_3 & p_2\\ p_2 & p_3 & p_0 & -p_1\\ p_3 & -p_2 & p_1 & p_0\\ \end{bmatrix} \begin{bmatrix} q_0 \\ q_1 \\ q_2 \\ q_3 \\ \end{bmatrix} \ \ \ (3)$

  • 「$r,p,q$」を以下のようなベクトル「$\vec{r},\vec{p},\vec{q}$」の成分からなる列ベクトルとする。
  • $\Large \displaystyle r= \begin{bmatrix} r_1 \\ r_2 \\ r_3 \\ \end{bmatrix} ,\ \ p= \begin{bmatrix} p_1 \\ p_2 \\ p_3 \\ \end{bmatrix} ,\ \ q= \begin{bmatrix} q_1 \\ q_2 \\ q_3 \\ \end{bmatrix} $

  • 「$\tilde{p}$」を以下のような外積行列とする。
  • $\Large \displaystyle \tilde{p}= \begin{bmatrix} 0 & -p_3 & p_2\\ p_3 & 0 & -p_1 \\ -p_2 & p_1 & 0 \\ \end{bmatrix} $

  • (3)式を「$r,p,q,\tilde{p}$」を使って表すと以下のようになる。
  • $\Large \displaystyle \begin{bmatrix} r_0 \\ r \\ \end{bmatrix} = \begin{bmatrix} p_0 & -p^T\\ p & p_0 I+\tilde{p} \\ \end{bmatrix} \begin{bmatrix} q_0\\ q \\ \end{bmatrix} $

4.2 Quaternionの性質

交換則
  • クオータニオンには、交換則が成立しない。
  • 例えば、次の演算結果は、
  • $\Large \displaystyle p*q - q*p= \begin{bmatrix} p_0 q_0-\vec{p}\cdot \vec{q}\\ p_0\vec{q}+q_0 \vec{p}+\vec{p}\times \vec{q} \\ \end{bmatrix} - \begin{bmatrix} q_0 p_0-\vec{q}\cdot \vec{p}\\ q_0\vec{p}+p_0 \vec{q}+\vec{q}\times \vec{p} \\ \end{bmatrix} = \begin{bmatrix} 0\\ 2\vec{q}\times \vec{p} \\ \end{bmatrix} $

  • となり、クオータニオンの積の順序は、交換できない。
  • $\Large \displaystyle p*q \neq q*p$

ノルム
積のノルム
共役クオータニオン
クオータニオンの逆数
積の逆数

4.3 Quaternionを使用したベクトルの回転

4.4 Quaternionの微分

クオータニオンの時間微分
  • クオータニオンの時間微分は、次の式で計算されます。
  • $\Large \displaystyle \dot{q}=\frac{1}{2}\omega q$

4.5 応用が広がるきっかけになった公式

応用が広がるきっかけになった公式
  • ある単位ベクトル「$u=(u_x, u_y, u_z)$」を軸として回転角「$\theta$」だけ右ねじの方向に回転させるとする.
  • この時,点「$A(x,y,z)$」の回転後の座標「$A'(x',y',z')$」はどうやって求めたらいいだろうか?それが四元数を使えば簡単に表せるというのである.
  • まず,次のような四元数を用意してやる.
  • $\Large \displaystyle X=\cos\frac{\theta}{2} + \sin\frac{\theta}{2}(u_x i + u_y j + u_z k)$

  • これが回転軸の向きと回転角度を表す四元数なわけだ.そして,点「$A$」の座標を表す四元数を次のように用意してやる.
  • $\Large \displaystyle A=x i + y j + z k$

  • 実数部分は 0 にしておく.そして次のような計算をしてやる.
  • $\Large \displaystyle B=X A \bar{X}$

  • 結果として得られる四元数「$B$」の実部は「$0$」になり,三つの虚部の方にはそれぞれ「$(x',y',z')$」の座標値が入っているというのである.これは驚きだ.

5.vector_mathライブラリの概要

5.1 vector_math

Aabb2
Defines a 2-dimensional axis-aligned bounding box between a min and a max position.
Aabb3
Defines a 3-dimensional axis-aligned bounding box between a min and a max position.
Colors
Contains functions for converting between different color models and manipulating colors. In addition to that, some known colors can be accessed for fast prototyping.
Frustum
Defines a frustum constructed out of six Planes.
IntersectionResult
Defines a result of an intersection test.
Matrix2
2D Matrix. Values are stored in column major order.
Matrix3
3D Matrix. Values are stored in column major order.
Matrix4
4D Matrix. Values are stored in column major order.
Obb3
Defines a 3-dimensional oriented bounding box defined with a center, halfExtents and axes.
Plane
Quad
Defines a quad by four points.
Quaternion
Defines a Quaternion (a four-dimensional vector) for efficient rotation calculations.
Ray
Defines a Ray by an origin and a direction.
Sphere
Defines a sphere with a center and a radius.
Triangle
Defines a triangle by three points.
Vector
Base class for vectors
Vector2
2D column vector.
Vector3
3D column vector.
Vector4
4D column vector.

5.2 vector_math_64

Aabb2
Defines a 2-dimensional axis-aligned bounding box between a min and a max position.
Aabb3
Defines a 3-dimensional axis-aligned bounding box between a min and a max position.
Colors
Contains functions for converting between different color models and manipulating colors. In addition to that, some known colors can be accessed for fast prototyping.
Frustum
Defines a frustum constructed out of six Planes.
IntersectionResult
Defines a result of an intersection test.
Matrix2
2D Matrix. Values are stored in column major order.
Matrix3
3D Matrix. Values are stored in column major order.
Matrix4
4D Matrix. Values are stored in column major order.
Obb3
Defines a 3-dimensional oriented bounding box defined with a center, halfExtents and axes.
Plane
Quad
Defines a quad by four points.
Quaternion
Defines a Quaternion (a four-dimensional vector) for efficient rotation calculations.
Ray
Defines a Ray by an origin and a direction.
Sphere
Defines a sphere with a center and a radius.
Triangle
Defines a triangle by three points.
Vector
Base class for vectors
Vector2
2D column vector.
Vector3
3D column vector.
Vector4
4D column vector.

5.3 vector_math_geometry

BarycentricFilter
CircleGenerator
ColorFilter
CubeGenerator
CylinderGenerator
FlatShadeFilter
GeometryFilter
GeometryGenerator
GeometryGeneratorFlags
InplaceGeometryFilter
InvertFilter
MeshGeometry
RingGenerator
SphereGenerator
TransformFilter
VertexAttrib

5.4 vector_math_lists

5.5 vector_math_operation

6.Quaternion classの概要

目次

6.1 Constructors

Quaternion(double x, double y, double z, double w)
  • Constructs a quaternion using the raw values x, y, z, and w. factory

Quaternion.axisAngle(Vector3 axis, double angle)
Constructs a quaternion from a rotation of angle around axis. factory

Quaternion.copy(Quaternion original)
Constructs a quaternion as a copy of original. factory

Quaternion.dq(Quaternion q, Vector3 omega)
Constructs a quaternion from time derivative of q with angular velocity omega. factory

Quaternion.euler(double yaw, double pitch, double roll)
  • Constructs a quaternion from yaw, pitch and roll.factory
  • 「オイラー角(yaw, pitch, roll)」から「Quaternion」を生成する。

Quaternion.fromBuffer(ByteBuffer buffer, int offset)
Constructs a quaternion with a storage that views given buffer starting at offset. offset has to be multiple of Float64List.bytesPerElement.

Quaternion.fromFloat64List(Float64List _qStorage)
Constructs a quaternion with given Float64List as storage.

Quaternion.fromRotation(Matrix3 rotationMatrix)
Constructs a quaternion from a rotation matrix rotationMatrix. factory

Quaternion.fromTwoVectors(Vector3 a, Vector3 b)
Constructs a quaternion to be the rotation that rotates vector a to b. factory

Quaternion.identity()
Constructs a quaternion set to the identity quaternion. factory

Quaternion.random(Random rn)
Constructs a quaternion with a random rotation. The random number generator rn is used to generate the random numbers for the rotation. factory

6.2 Properties

axis → Vector3
  • axis of rotation.
  • read-only

6.3 Methods

absoluteError(Quaternion correct) → double
  • Absolute error between this and correct.

add(Quaternion arg) → void
  • Add arg to this.

asRotationMatrix() → Matrix3
  • Returns a rotation matrix containing the same rotation as this.

6.4 Operators