クォータニオンによる回転の原理をスッキリ理解したい

※ 本記事は、もともと Qiita に書いた記事の転載です

クォータニオンによる回転のモヤモヤ

ここではクォータニオンの定義や演算規則については書きません。 いきなりクォータニオンによる3次元座標の回転変換について結果を眺めるところから話を始めましょう。

3次元座標 \((x, y, z)\) を次のようにクォータニオンで表します。(実部は何でも良いのですがここではゼロとします)

\[p = ix + jy + kz\]

また、回転軸 \(n\) の周りに角度 \(\theta\) 回転する変換を表すクォータニオンは次式で表されます。

\[\begin{split}q = \cos\frac{\theta}{2} + \boldsymbol{i}^T \boldsymbol{n}\sin\frac{\theta}{2}, \quad \boldsymbol{i} = \left(\begin{array}{c} i \\ j \\ k\end{array}\right)\end{split}\]

これを用いて回転後の座標は次式で計算できます。

\[p' = qpq^*\]

これについて初めて学んだときの私の感想は次のようなものでした。

  • 式を展開すれば回転変換になっていることは確かめられたが、式の与えられ方が天下り的に感じられモヤモヤする。
  • 2次元平面の回転変換が複素数で表されることと関係がありそうな感じがするが、その関係がハッキリとは分からない。
  • 複素数のときは \(\theta\) の回転は \(e^{i\theta}\) で表されたのに、クォータニオンでは何故 \(\theta/2\) が出て来るのだろう?
  • 複素数のときは片側から掛ければよかったのに、クォータニオンでは何故両側から挟み撃ちするように掛ける必要があるのだろう?

当時、このモヤモヤをズバッと解消してくれる説明を私は見つけることが出来ませんでした。そこで試行錯誤して式をこねくり回し、自分なりの結論を得ました。この記事はその試行錯誤をまとめ直したものです。

Z軸周りの回転

まずは \(p\) をZ軸回りに回転する変換について考えていきます。そのために \(p\) を次のように変形します。

\[\begin{split}\begin{array}{cl} p &= ix + jy + kz \\ &= (x + ky)i + kz \end{array}\end{split}\]

この式変形のポイントを箇条書しておきます。

  • 途中 \(ki = j\) というクォータニオンの演算ルールを用いて \(j\) を消しました。
  • 第1項 \(x + ky\)\(k\) を虚数単位とする複素数となっています。
  • 今はZ軸回りの回転を考えたいので、第2項を保ったまま第1項をXY平面上で回転できれば良さそうです。

さて、第1項は複素数なので、次のように指数関数表示に書き換えることができます。

\[x + ky = \sqrt{x^2 + y^2}(\cos\alpha + k\sin\alpha) = ae^{k\alpha}\]

これを代入すると \(p\) は次式となります。

\[p = \left(ae^{k\alpha}\right)i + kz\]
さあ、この形になると、なんだか複素数の知識を利用して回転変換ができそうな気がしてきませんか?してきましたね?
\(p\) をZ軸回りに \(\theta\) 回転した値は次式で表せます。
\[p' = \left(ae^{k(\alpha + \theta)}\right)i + kz\]

どういった演算をすれば \(p\) から \(p'\) を得ることが出来るでしょうか。 複素数からの雑な連想ですが、つぎのクォータニオンを掛けるとどうなるでしょうか。

\[q = \cos\theta + k\sin\theta = e^{k\theta}\]

では試しにやってみましょう。左から掛けてみます。

\[\begin{split}\begin{array}{cl} qp &= e^{k\theta}\left\{\left(ae^{k\alpha}\right)i + kz \right\} \\ &= \left(ae^{k(\alpha + \theta)}\right)i + e^{k\theta}kz \end{array}\end{split}\]

あれれ、失敗です。 まず第1項は目論見どおりで成功です。第1項は x, y 座標を複素数で表現している項なので、この偏角に \(\theta\) が加わっているのはXY平面上での回転を表しています。 しかし第2項の値が変化してしまっているのは頂けませんね。第2項はZ座標を表す項です。Z軸回りの回転では第2項は変化せずに保たれなければなりません。

これを解決するため、まず \(qp\) の式を次のように変形します。

\[qp = i\left(ae^{-k(\alpha + \theta)}\right) + e^{k\theta}kz\]

第1項の \(i\) を右から左へ移動し、偏角にマイナス符号が加わっています。この変形は \(ki = -ik\) というクォータニオンの演算ルールから導けます。 これに対して右から \(q^* = e^{-k\theta}\) を掛けてみます。

\[\begin{split}\begin{array}{cl} qpq^* &= \left\{i\left(ae^{-k(\alpha + \theta)}\right) + e^{k\theta}kz \right\}e^{-k\theta} \\ &= i\left(ae^{-k(\alpha + \theta) - k\theta}\right) + e^{k\theta}kze^{-k\theta} \\ &= i\left(ae^{-k(\alpha + 2\theta)}\right) + kz \end{array}\end{split}\]

Z軸回りに \(2\theta\) 回転した式になりました! 結局、q

\[q = \cos\frac{\theta}{2} + k\sin\frac{\theta}{2} = e^{k\frac{\theta}{2}}\]

と置き直せば \(qpq^*\) がZ軸回りの \(\theta\) 回転を表すことが分かります。これでかなりのモヤモヤが晴れました。

  • Z軸回りの回転を、複素数の回転の知識を利用しつつ導出することが出来た。
  • 片側から単純に掛けるとZ成分(第2項)も変化してしまうため、これをキャンセルするために両側から \(\theta/2\) ずつ挟み撃ちする必要があることが分かった。

ただし、今までの議論はZ軸回りの回転に限定していました。この議論を任意の軸回りの回転に拡張していきましょう。

任意の軸回りの回転

任意の回転軸ベクトル \(n\) (単位ベクトル)の周りの回転を考えていきます。 次の順序で進めます。

  1. \(n\) をZ軸とする正規直交基底 \(u, v, n\) を作ります。
  2. 虚数単位 \(i, j, k\)\(\{u, v, n\}\)-座標系に座標変換します。
  3. その座標系でZ軸回りに回転すれば \(n\) 軸回りの回転が得られます。

(1) 正規直交基底を作る

与えられた回転軸 \(n\) を使って、正規直交基底 \(u, v, n\) を構成します。もちろん \(u, v\) は一意には定まりませんが、正規直交基底の条件を満たしていればなんでも構いません。 正規直交基底なので、u, v, n には次の関係が成り立ちます。

\[\begin{split}||u|| = ||v|| = ||n|| = 1 \\ u\cdot v = v\cdot n = n\cdot u = 0 \\ u \times v = n, \quad v \times n = u, \quad n \times u = v\end{split}\]

これらの基底ベクトルを並べた直交行列 \(U\) を定義しておきます。

\[U=(\array{u & v & n})\]

これが座標変換行列となります。

(2) 虚数単位を座標変換する

\[\begin{split}\boldsymbol{i} = \left(\array{i \\ j \\ k}\right), \quad \boldsymbol{r} = \left(\array{x \\ y \\ z}\right)\end{split}\]

と置くと、 \(p\) は次のように書くことが出来ます。

\[p = ix + jy + kz = \boldsymbol{i}^T\boldsymbol{r}\]

さて、前項の正規直交座標系への座標変換は次式で表されます。

\[r = Ur'\]

これを代入して式変形を行います。

\[\begin{split}\begin{array}{ll} p &= i^Tr \\ &= i^TUr' \\ &= (\array{i^Tu & i^Tv & i^Tn})r' \\ &= (\array{i' & j' & k'})r' \\ &= i'x' + j'y' + k'z' \end{array}\end{split}\]

ただし、

\[\begin{split}i' = \left(\array{i' \\ j' \\ k'}\right) = \left(\array{i^Tu \\ i^Tv \\ i^Tn}\right)\end{split}\]

なんと実は \(i', j', k'\) も下記のクォータニオン虚数単位の性質を満たすことが示されます(導出は後述)

\[\begin{split}i'^2 = j'^2 = k'^2 = -1 \\ i'j' = -j'i' = k', \quad j'k' = -k'j' = i', \quad k'i' = -i'k' = j'\end{split}\]

つまり、座標変換行列 \(U\) によって虚数単位 \(i, j, k\) を別の虚数単位 \(i', j', k'\) に座標変換することが出来ました!

(3) Z’軸回りに回転

ここまでくればもう簡単です。ベクトル \(n\) 周りに回転するためには、新しい座標系で Z’ 軸回りに回転すれば良いのです。

\(p\) は新しい虚数単位 \(i', j', k'\) を使って次の形に式変形できます。
\[p = i'x' + j'y' + k'z' = \left(ae^{k\alpha}\right)i' + k'z'\]

これを Z’ 軸回りに回転するクォータニオンは次式で与えられるのでした。

\[q = \cos\frac{\theta}{2} + k'\sin\frac{\theta}{2}\]

ここで

\[k' = i^Tn\]

でしたので、これを代入することにより次式が得られます。

\[q = \cos\frac{\theta}{2} + i^Tn\sin\frac{\theta}{2}\]

やった!目的の式が得られました。

虚数単位の座標変換の詳細

落ち穂拾いです。 ここでは座標変換で得られた新しい虚数単位 \(i', j', k'\) もクォータニオンの性質を満たすことを導出したいと思います。

準備として、次の2つのクォータニオンの掛け算について調べていきます。

\[\begin{split}q_1 = i x_1 + jy_1 + kz_1 = i^Tv_1, \quad v_1 = (\array{x_1 & y_1 & z_1})^T \\ q_2 = i x_2 + jy_2 + kz_2 = i^Tv_2, \quad v_2 = (\array{x_2 & y_2 & z_2})^T \\\end{split}\]

これらの掛け算は次式となります。

\[q_1 q_2 = i^Tv_1\ i^Tv_2 = v_1^T(ii^T)v_2\]

間に挟まれている行列 \(ii^T\) が重要です。

\[\begin{split}ii^T = \left(\array{i \\ j \\ k}\right)(\array{ i & j & k }) = \left(\array{ -1 & k & -j \\ -k & -1 & i \\ j & -i & -1 }\right)\end{split}\]

ここで \(I\) を単位行列とし、次の行列 \(K\) を定義すれば

\[\begin{split}K= \left(\array{ 0 & k & -j \\ -k & 0 & i \\ j & -i & 0 }\right)\end{split}\]

\(ii^T\) は次のように表せます。

\[ii^T = K - I\]

これを代入すると先程のクォータニオンの積は次式となります。

\[q_1 q_2 = v_1^T(K - I)v_2 = v_1^TKv_2 - v_1\cdot v_2\]

さてここで、これは式を展開すれば分かるのですが、次式のように第1項からは外積ベクトルが得られます。

\[v_1^TKv_2 = i^T(v_1\times v_2)\]

以上で準備は終わりです。いよいよ座標変換された虚数単位 \(i', j', k'\) の性質を調べていきましょう。

\[i'^2 = (i^Tu)^2 = u^T(K-I)u = i^T(u \times u) - u\cdot u = -1\]
\[i'j' = (i^Tu)(i^Tv) = u^T(K-I)v = i^T(u \times v) - u\cdot v = i^Tn = k'\]

以下同様にして、クォータニオンの虚数単位の性質を確かめられます。