3次元回転を扱う、ロドリゲスの回転公式とクォータニオンの中間の方法

3次元回転を扱う際、以下のような方法がよく知られています。

1. オイラー角
2. 回転行列
3. ロドリゲスの回転公式（軸角表現）
4. クォータニオン

これらに比べて言及されることが少ないですが、「[オイラー–ロドリゲスの公式](https://en.wikipedia.org/wiki/Euler–Rodrigues_formula)」というものがあります。知っておくと便利だと感じたため纏めることにしました。

## 公式（ベクトル形式）

※行列形式については「[回転行列との相互変換](#回転行列との相互変換)」の節で触れます。

ロドリゲスの回転公式から簡単に変形できるため、やっていきます。

軸角表現では名前通り、回転軸を表す3次元単位ベクトル $\vec{k} = (x, y, z)$ と回転角 $\theta$ を用意します。これらを用いて、ロドリゲスの回転公式は以下の式で表されます。

```math
\begin{align}
\vec{v}' &= \vec{v} \cos{\theta} + \left( \vec{k} \times \vec{v} \right) \sin{\theta} + \vec{k} \left( \vec{k} \cdot \vec{v} \right) \left( 1 - \cos{\theta} \right) \\
&= \vec{v} + \left( \vec{k} \times \vec{v} \right) \sin{\theta} + \left( \vec{k} \times \left( \vec{k} \times \vec{v} \right) \right) \left( 1 - \cos{\theta} \right) \\
\end{align}
```

ここで $\theta = 2 \cdot \frac{\theta}{2}$ と置き換えて2倍角の公式を適用すると、 $\sin{\theta} = 2 \sin{\frac{\theta}{2}} \cos{\frac{\theta}{2}}$, $(1 - \cos{\theta}) = 2 \sin^2{\frac{\theta}{2}}$ となります。すると、ちょうど $\vec{k}$ と $\sin{\frac{\theta}{2}}$ の数が各項で一致するため、これらは纏めることができます。

$a = \cos{\frac{\theta}{2}}$, $\vec{\omega} = (b, c, d) = \sin{\frac{\theta}{2}}$ として整理すると、以下の「オイラー–ロドリゲスの公式」が得られます。三角関数が式から消えました。

```math
\vec{v}' = \vec{v} + 2a \left( \vec{\omega} \times \vec{v} \right) + 2 \left( \vec{\omega} \times \left( \vec{\omega} \times \vec{v} \right) \right)
```

4つのパラメータを改めて書くと以下の通りです（「オイラーパラメータ」などと呼ぶそうです）。 $\theta, x, y, z$ は軸角表現のものであるため、回転軸と回転角からパラメータをすぐ作成できることがわかります。

```math
(a, b, c, d) = \left( \cos{\frac{\theta}{2}}
, \; x \sin{\frac{\theta}{2}}
, \; y \sin{\frac{\theta}{2}}
, \; z \sin{\frac{\theta}{2}} \right)
```

パラメータの主な特徴は以下のとおりです。

* 4成分の二乗和は 1
* そのうち3成分は回転軸に平行
* 無回転で回転軸を決められない状況でも問題なく表現できる
* 欠点？は「 $\theta$ が 360° 異なる回転」を考えるとわかるとおり、同じ回転を表すパラメータが正負の2通りある（補間の際に罠となる）

---

さて、パラメータの定義を見れば気付くと思いますが、パラメータは**回転を表す単位クォータニオンの成分と全く同じ**です。この公式を用いれば、クォータニオンで回転するあの式を真面目に計算する必要は無く、したがってクォータニオンの四則演算などを実装する必要もありません。

### クォータニオンによる回転の式からの導出

パラメータが単位クォータニオンの成分と同じということは、クォータニオンによる回転の式 $\vec{v}' = q \vec{v} \bar{q}$ を変形しても同じ公式を導けるはずです。それを軽く確認しておきます。

<details><summary>（クリックして展開）</summary><div>
簡単に計算するためクォータニオンは「スカラー＋3次元ベクトル」で表現し、乗算記号が無いベクトル同士の積は $\vec{v_1} \vec{v_2} = - \vec{v_1} \cdot \vec{v_2} + \vec{v_1} \times \vec{v_2}$ と変形していきます。回転を表す単位クォータニオンは $q = a + \vec{\omega}$ とします。

```math
\begin{align}
\vec{v}' &= q \vec{v} \bar{q} \\
&= \left( a + \vec{\omega} \right) \left( 0 + \vec{v} \right) \left( a - \vec{\omega} \right) \\
&= \left( a \vec{v} - \vec{\omega} \cdot \vec{v} + \vec{\omega} \times \vec{v} \right) \left( a - \vec{\omega} \right) \\
&= a^2 \vec{v} - \cancel{a \vec{\omega} \cdot \vec{v}} + a \vec{\omega} \times \vec{v} \\
&\quad\; + \cancel{a \vec{v} \cdot \vec{\omega}} - a \vec{v} \times \vec{\omega} + \left( \vec{\omega} \cdot \vec{v} \right) \vec{\omega} + \cancel{\left( \vec{\omega} \times \vec{v} \right) \cdot \vec{\omega}} - \left( \vec{\omega} \times \vec{v} \right) \times \vec{\omega} \\
&= \left( 1 - \vec{\omega} \cdot \vec{\omega} \right) \vec{v} + 2 a \left( \vec{\omega} \times \vec{v} \right) + \vec{\omega} \left( \vec{\omega} \cdot \vec{v} \right) + \left( \vec{\omega} \times \left( \vec{\omega} \times \vec{v} \right) \right) \\
&= \vec{v} + 2 a \left( \vec{\omega} \times \vec{v} \right) + 2 \left( \vec{\omega} \times \left( \vec{\omega} \times \vec{v} \right) \right) \\
\end{align}
```

きちんとスカラー部が消え、オイラー–ロドリゲスの公式の形に整理できました。
</div></details>

## 応用

パラメータがクォータニオンのものと同じであるということは、クォータニオンでできることはこちらの方法にも適用できます。

（というより逆に考えて、普段は回転の情報をクォータニオンで管理・計算しておき、最後に回転を適用するときだけ成分を取り出して使うとみたほうが自然かもしれません）

### 回転の合成

クォータニオンの積と同じです。（具体的な式は省略します）

任意の回転の合成は、オイラー角や軸角表現では難しく、回転行列かクォータニオンが適しています。特にクォータニオンであれば4成分で済むため簡潔です。

### 回転行列との相互変換

ベクトルの式を回転行列 $R$ に変形すると以下のようになります。こちらも三角関数が現れません。

```math
\vec{v}' = R \vec{v} = \begin{pmatrix}
a^2 + b^2 - c^2 - d^2 & 2(bc - ad) & 2(bd + ac) \\
2(bc + ad) & a^2 + c^2 - b^2 - d^2 & 2(cd - ab) \\
2(bd - ac) & 2(cd + ab) & a^2 + d^2 - b^2 - c^2 \\
\end{pmatrix} \vec{v}
```

逆に回転行列からパラメータを求めるのは、例えば以下の式で実行できます。

```math
\begin{align}
a &= \frac{1}{2} \sqrt{1 + R_{11}^2 + R_{22}^2 + R_{33}^2} \\
b &= \frac{1}{4a} \left( R_{32} - R_{23} \right) \\
c &= \frac{1}{4a} \left( R_{13} - R_{31} \right) \\
d &= \frac{1}{4a} \left( R_{21} - R_{12} \right) \\
\end{align}
```

考え方としては、2乗和が 1 であることおよび行列の対角成分から、各パラメータの2乗がすぐ求まります。しかし2乗を外すときに符号が矛盾しないよう決める必要があるため、ひとつだけ符号を決めて、残りは非対角成分から計算します。
※それでも符号の決め方が2通りとなりますが、これは同じ回転を表すパラメータが正負の2通りあるためです。

### 回転の補間

4成分であるパラメータは4次元単位ベクトルと見なせ、4次元空間中にプロットすれば3次元球面上の点と対応付けられます。球面上で2点の間の点を求めることで、中間の回転を表すパラメータを得られます。この補間は「[球面線形補間](https://en.wikipedia.org/wiki/Slerp)」（spherical linear interpolation; slerp）で均等に行えます。

具体的には、2つの単位ベクトル $\vec{p_1}, \vec{p_2}$ を $t:(1-t)$ で補間する場合、以下のようになります。

```math
\vec{p} = \frac{\sin{\left( (1-t) \Omega \right)}}{\sin{\Omega}} \vec{p_1} + \frac{\sin{\left( t \Omega \right)}}{\sin{\Omega}} \vec{p_2}
\quad \textrm{where} \quad \Omega = \arccos{\left( \vec{p_1} \cdot \vec{p_2} \right)}
```

均等でなくてもよいなら、単純な線形補間で求めたベクトルを長さ 1 となるよう規格化する方法もあります。

なお、符号反転した $-\vec{p_2}$ も同じ回転を表すため、 $\vec{p_1}$ に近いほうをゴールに選んで補間しないと奇妙な回転になってしまいます。近いかどうかは $\Omega$ を求める際の内積の符号で判定できます。

## まとめ

3次元回転を扱う「オイラー–ロドリゲスの公式」について整理しました。

これは**クォータニオンの方法と同じパラメータ**で計算しますが、クォータニオンのまま計算するわけではないため、クォータニオンを実装する必要はありません。とはいえパラメータが同じため、回転の合成や球面線形補間などクォータニオンが有利な部分を真似ることができます。

ロドリゲスの回転公式と比較すると、三角関数がパラメータに内包されていて公式側には現れないという特徴があります。回転軸や回転角を知る必要が無いときは、三角関数を用いる手間を省けます。