Quaternion

앞선 포스팅에서 살펴보았듯, Rotation matrix는 redundancy를 가지고 있고 Rotation vector나 Euler angle은 singularity를 가지고 있다. 이러한 문제를 해결하기 위해, 단위 복소수를 이용하여 2D 평면에서 회전을 표현했던 것을 이용할 수 있다. 아래 수식과 같이 $\theta$만큼의 회전을 복소수로 표현할 수 있다.

$e^{i\theta} = \cos\theta + i\sin\theta$

3차원 공간에서는 사원수라고 불리는 Quaternion을 이용하여 회전을 표현할 수 있다. Quaternion은 실수부가 1개, 허수부가 3개인 수이다. 아래 수식에서 $i$, $j$, $k$는 허수부를 나타낸다. 모든 허수부가 0이면 real quaternion, 실수부가 0이면 imaginary quaternion이라고 부른다.

$\mathbf{q} = s + xi + yj + zk$

Quaternion은 아래와 같이 크기가 4인 벡터로도 표현할 수 있다.

$\mathbf{q} = [s, \mathbf{v}]^T, \quad s = s \in \mathbb{R}, \quad \mathbf{v} = [x, y, z]^T \in \mathbb{R}^3$

Quaternion의 세 허수부는 아래 수식과 같은 규칙을 따른다.

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

제곱을 할 경우 복소수와 마찬가지로 -1이 되며, $i$, $j$, $k$를 각각 $x$, $y$, $z$축으로 본다면 각 허수부끼리의 곱셈은 각 축에 해당하는 방향벡터끼리의 외적으로 볼 수 있다.

Quaternion Operations

이 파트에서는 사원수간 연산 규칙에 대해 다루도록 하겠다. 일반적인 복소수와 비슷한 규칙을 가진다.

Addition and Subtraction

$\mathbf{q}_1 \pm \mathbf{q}_2 = [s_1 \pm s_2, \mathbf{v}_1 \pm \mathbf{v}_2]^T$

Multiplication

두 사원수 $\mathbf{q}_1 = [s_1, \mathbf{v}_1]^T$와 $\mathbf{q}_2 = [s_2, \mathbf{v}_2]^T$의 곱셈은 일반 복소수와 같이 분배법칙으로 계산할 수 있다.

$\begin{align*} \mathbf{q}_1\mathbf{q}_2 & = s_1s_2 - x_1x_2 - y_1y_2 - z_1z_2 \\ & + (s_1x_2 + x_1s_2 + y_1z_2 - z_1y_2)i \\ & + (s_1y_2 - x_1z_2 + y_1s_2 + z_1x_2)j \\ & + (s_1z_2 + x_1y_2 - y_1x_2 + z_1s_2)k \end{align*}$

아래와 같이 벡터 연산 형태로 간단하게 표현할 수도 있다.

$\mathbf{q}_1\mathbf{q}_2 = [s_1s_2 - \mathbf{v}_1^T\mathbf{v}_2, s_1\mathbf{v}_2 + s_2\mathbf{v}_1 + \mathbf{v}_1 \times \mathbf{v}_2]^T$

이때, 벡터 연산 수식을 보면 외적이 포함되어 있기 때문에 $\mathbf{v}_1$, $\mathbf{v}_2$가 평행하지 않다면 교환법칙은 성립하지 않는다. 즉,

$\mathbf{q}_1\mathbf{q}_2 \neq \mathbf{q}_2\mathbf{q}_1$.

Quaternion 곱셈의 결합 법칙은 성립한다.

$(\mathbf{q}_1\mathbf{q}_2)\mathbf{q}_3 = \mathbf{q}_1(\mathbf{q}_2\mathbf{q}_3)$

Length

Quaternion의 크기는 다음과 같이 정의된다.

$\Vert\mathbf{q}\Vert = \sqrt{s^2 + x^2 + y^2 + z^2}$

이때 두 quaternion의 곱의 크기는 각 quaternion의 크기를 곱한 것과 같다. 이 성질 덕분에 unit quaternion을 계속 곱하더라도 그 크기가 유지된다.

$\Vert\mathbf{q}_1\mathbf{q}_2\Vert = \Vert\mathbf{q}_1\Vert\Vert\mathbf{q}_2\Vert$

Conjugate

Quaternion의 켤레는 다음과 같이 모든 허수부의 부호를 바꾸는 것으로 정의된다.

$\mathbf{q}^* = s - xi - yj - zk = [s, -\mathbf{v}]^T$

어떤 복소수와 그 켤레를 곱했을 경우 실수부가 그 복소수의 크기의 제곱인 real quaternion이 된다.

$\mathbf{q}\mathbf{q}^* = \mathbf{q}^*\mathbf{q} = \Vert\mathbf{q}\Vert^2$

Inverse

Quaternion의 역수은 곱해서 1이 되는 수이다.

$\mathbf{q}\mathbf{q}^{-1} = \mathbf{q}^{-1}\mathbf{q} = 1$

역수는 다음과 같이 계산할 수 있다.

$\mathbf{q}^{-1} = \frac{\mathbf{q}^*}{\Vert\mathbf{q}\Vert^2}$

이때, unit quaternion의 경우 켤레가 역수와 같기 때문에 unit quaternion을 곱한 뒤 역수를 취하면 행렬의 역을 취한 것과 유사한 결과가 나온다.

$(\mathbf{q_1}\mathbf{q_2})^{-1} = \mathbf{q_2}^{-1}\mathbf{q_1}^{-1}$

Scalar Multiplication

Quaternion의 상수배는 일반 복소수와 같이 정의된다.

$c\mathbf{q} = [cs, c\mathbf{v}]^T$

Use Quaternion to Represent a Rotation

Unit quaternion을 이용하여 어떻게 3차원 회전을 계산할 수 있을지 살펴보자.먼저 3차원 공간상의 점 $\mathbf{p} = [x, y, z]^T$가 있다고 하자. 우선 크기를 맞추기 위해 좌표값 3개를 허수부에 넣고 실수부에는 0을 할당하여 $\mathbf{p} = [0, x, y, z]^T$로 변환한다. 회전된 점을 $\mathbf{p}'$라고 하면, unit quaternion을 이용하여 아래와 같이 연산을 수행할 수 있다.

$\mathbf{p}' = \mathbf{q}\mathbf{p}\mathbf{q}^{-1} = \mathbf{q}\mathbf{p}\mathbf{q}^*$

이때의 곱셈은 위에서 정의한 quaternion간의 곱셈이고, 그 결과인 $\mathbf{p}'$는 실수부가 0인 unit quaternion이다. 따라서 $\mathbf{p}'$의 허수부를 취하여 변환된 3차원 점을 구할 수 있다.

Conversion of Quaternions to Other Rotation Representations

Quaternion을 rotation matrix와 rotation vector로 변환하는 방법에 대해 살펴보자. 우선 두 가지 연산자 $^+$와 $^\oplus$를 정의하자.

$\mathbf{q} = [s, x, y, z]^T = [s, \mathbf{v}]^T$

일 때,

$\mathbf{q}^+ = \begin{bmatrix} s & -\mathbf{v}^T \\ \mathbf{v} & s\mathbf{I} + \mathbf{v}^\wedge \end{bmatrix} = \begin{bmatrix} s & -x & -y & -z \\ x & s & -z & y \\ y & z & s & -x \\ z & -y & x & s \end{bmatrix}$

$\mathbf{q}^\oplus = \begin{bmatrix} s & -\mathbf{v}^T \\ \mathbf{v} & s\mathbf{I} - \mathbf{v}^\wedge \end{bmatrix} = \begin{bmatrix} s & -x & -y & -z \\ x & s & z & -y \\ y & -z & s & x \\ z & y & -x & s \end{bmatrix}$

위 두 연산은 quaternion을 4x4 행렬로 변환시킨다. 두 연산자 $^+$와 $^\oplus$를 사용하면 quaternion 곱셈을 아래와 같이 표현할 수 있다. 이때 $\mathbf{q}_1\mathbf{q}_2$는 quaternion 간 곱셈, 나머지 둘은 행렬 간 곱셈이다.

$\mathbf{q}_1\mathbf{q}_2 = \mathbf{q}_1^+ \mathbf{q}_2 = \mathbf{q}_2^\oplus \mathbf{q}_1$

각 연산자를 풀어서 쓰면 아래와 같다.

$\mathbf{q}_1^+\mathbf{q}_2 = \begin{bmatrix} s_1 & -\mathbf{v}_1^T \\ \mathbf{v}_1 & s_1\mathbf{I} + \mathbf{v}_1^\wedge \end{bmatrix} \begin{bmatrix} s_2 \\ \mathbf{v}_2 \end{bmatrix} = \begin{bmatrix} s_1s_2 - \mathbf{v}_1^T\mathbf{v}_2 \\ s_1\mathbf{v}_2 + s_2\mathbf{v}_1 + \mathbf{v}_1^\wedge\mathbf{v}_2 \end{bmatrix} = \mathbf{q}_1\mathbf{q}_2$

$\mathbf{q}_2^\oplus\mathbf{q}_1 = \begin{bmatrix} s_2 & -\mathbf{v}_2^T \\ \mathbf{v}_2 & s_2\mathbf{I} - \mathbf{v}_2^\wedge \end{bmatrix} \begin{bmatrix} s_1 \\ \mathbf{v}_1 \end{bmatrix} = \begin{bmatrix} s_2s_1 - \mathbf{v}_2^T\mathbf{v}_1 \\ s_2\mathbf{v}_1 + s_1\mathbf{v}_2 - \mathbf{v}_2^\wedge\mathbf{v}_1 \end{bmatrix} = \mathbf{q}_1\mathbf{q}_2$

Quaternion to Rotation Matrix

연산자를 정의했으니 이제 quaternion을 이용한 회전변환 식을 rotation matrix를 이용한 식으로 바꾸어 보자.

$\mathbf{p}' = \mathbf{q}\mathbf{p}\mathbf{q}^{-1}$

위 식을 $^+$ 연산자를 이용하여 행렬 연산으로 바꿀 수 있다.

$\mathbf{p}' = \mathbf{q}^+\mathbf{p}^+\mathbf{q}^{-1}$

이후 $\mathbf{q}_1\mathbf{q}_2 = \mathbf{q}_2^\oplus \mathbf{q}_1$를 이용하여 $p^+$와 $q^{-1}$의 위치를 바꾼다.

$\mathbf{p}' = \mathbf{q}^+{\mathbf{q}^{-1}}^\oplus\mathbf{p}$

$\mathbf{q}^+{\mathbf{q}^{-1}}^\oplus$ 부분을 풀어서 정리하면 아래와 같다.

$\mathbf{q}^+{\mathbf{q}^{-1}}^\oplus = \begin{bmatrix} s & -\mathbf{v}^T \\ \mathbf{v} & s\mathbf{I} + \mathbf{v}^\wedge \end{bmatrix} \begin{bmatrix} s & \mathbf{v}^T \\ -\mathbf{v} & s\mathbf{I} + \mathbf{v}^\wedge \end{bmatrix} = \begin{bmatrix} 1 & \mathbf{0} \\ \mathbf{0}^T & \mathbf{v}\mathbf{v}^T + s^2\mathbf{I} + 2s\mathbf{v}^\wedge + (\mathbf{v}^\wedge)^2 \end{bmatrix}$

이를 원래 수식에 대입하면 아래와 같다.

$\mathbf{p}' = \begin{bmatrix} 0 \\ x' \\ y' \\ z' \end{bmatrix} = \begin{bmatrix} 1 & \mathbf{0} \\ \mathbf{0}^T & \mathbf{v}\mathbf{v}^T + s^2\mathbf{I} + 2s\mathbf{v}^\wedge + (\mathbf{v}^\wedge)^2 \end{bmatrix} \begin{bmatrix} 0 \\ x \\ y \\ z \end{bmatrix}$

이때 첫 번째 row 또는 column은 계산에 영향을 미치지 못하므로 우하단 3개 열이나 행만을 보면 아래와 같다.

$\begin{bmatrix} x' \\ y' \\ z' \end{bmatrix} = \begin{bmatrix} \mathbf{v}\mathbf{v}^T + s^2\mathbf{I} + 2s\mathbf{v}^\wedge + (\mathbf{v}^\wedge)^2 \end{bmatrix} \begin{bmatrix} x \\ y \\ z \end{bmatrix}$

즉, 아래 공식을 통해 rotation matrix $\mathbf{R}$를 구할 수 있다.

$\mathbf{R} = \begin{bmatrix} \mathbf{v}\mathbf{v}^T + s^2\mathbf{I} + 2s\mathbf{v}^\wedge + (\mathbf{v}^\wedge)^2 \end{bmatrix}$

Quaternion to Rotation Vector

Quaternion으로부터 rotation vector를 구하기 위해서는 먼저 rotation matrix를 구한 뒤 trace를 취해야 한다. 위에서 구한 $\mathbf{R}$을 이용하여 수식을 정리하면 아래와 같다.

$\text{tr}(\mathbf{R}) = 4s^2 - 1$

또한 지난 포스팅에서는 rotation matrix를 rotation vector로 변환하는 과정에서 아래의 수식을 사용하였다.

$\theta = \arccos\left(\frac{\text{tr}(\mathbf{R}) - 1}{2}\right)$

위의 식을 아래 식에 대입하면

$\theta = \arccos\left(\frac{4s^2 - 1 - 1}{2}\right) = \arccos(2s^2 - 1)$

이때 arccos을 좌변으로 넘기고 반각공식을 사용하여 간단하게 정리할 수 있다.

$\cos\theta = 2s^2 - 1 = 2\cos^2\frac{\theta}{2} - 1$

$s = \cos\frac{\theta}{2}$

최종적으로 아래의 수식을 얻는다.

$\theta = 2\arccos(s)$

방향을 나타내는 unit vector $\mathbf{n}$은 quaternion $\mathbf{q}$의 허수부만을 뽑아 normalizing하는 것으로 구할 수 있다. 그 이유는, $\mathbf{q}$의 허수부를 3차원 점으로 삼아 $\mathbf{q}$로 회전시키면 같은 점이 나오는데, 그 말은 그 점이 회전축 위에 있다는 뜻이기 때문이다.

$\mathbf{n} = \frac{\mathbf{v}}{\Vert\mathbf{v}\Vert} = [q_1, q_2, q_3]^T / \sin\frac{\theta}{2}$