Quaternion
Quaternions are a number system which can be used to represent rotations in 3D space.
An extension called Dual quaternion, involving quaternions of dual numbers or two quaternions of real numbers, are able to represent both rotations and translations in 3D space.
Background
Quaternions are a number space that can be represented as a 4D point \(\displaystyle q = (q_0, q_1, q_2, q_3) = (q_0, \mathbf{q})\) with unit norm.
Here, \(\displaystyle \mathbf{q}\) represents the imaginary part: \(\displaystyle q_1 \mathbf{i}+q_2 \mathbf{j}+q_3 \mathbf{k}\).
The conjugate is \(\displaystyle \bar{q} = (q_0, -\mathbf{q})\).
Multiplication
Multiplication is \(\displaystyle q * p = (q_0 p_0 - \mathbf{q} \cdot \mathbf{p}, q_0 \mathbf{p} + p_0 \mathbf{q} + \mathbf{q} \times \mathbf{p})\).
Similar to imaginary numbers, \(\displaystyle \mathbf{i}^2=\mathbf{j}^2=\mathbf{k}^2=\mathbf{i}\mathbf{j}\mathbf{k}=-1\).
Multiplication between bases follow cross product rules: \(\displaystyle \mathbf{i} * \mathbf{j} = \mathbf{k}\) and \(\displaystyle \mathbf{j} * \mathbf{i} = -\mathbf{k}\).
In general, quaternion multiplication does not commute: \(\displaystyle q_1 * q_2 \neq q_2 * q_1\)
Intuition: Left multiplication by a quaternion is a 4D rotation plus a scaling operation.
Rotations
Quadratically conjugating a quaternion as follows is equivalent to applying a rotation on \(\displaystyle \mathbf{x}=(x,y,z)\):
\[\mathbf{q} * (c,x,y,z)*\bar{\mathbf{q}} = (c, R(\mathbf{q}) \cdot \mathbf{x})\]
Here, \(\displaystyle R(q)=R(-q)\) is a two-to-one mapping from quaternions to 3x3 rotation matrices.
Quaternion multiplication corresponds to composing rotations:
\(\displaystyle R(q*p) = R(q) \cdot R(p)\)
The quaternion \(\displaystyle q=(\cos(\theta/2), \hat{n}_1 \sin(\theta/2), \hat{n}_2 \sin(\theta/2), \hat{n}_3 \sin(\theta/2))\) is equivalent to the rotation around axis \(\displaystyle \hat{n}\) by angle \(\displaystyle \theta\).
Slurp
Slerp, or spherical linear interpolation is done as: \[ \operatorname{slerp}(q_0, q_1, s) \equiv q(s)[q_0, q_1] = q_0 \frac{\sin((1-s)\phi)}{\sin \phi} + q_1 \frac{\sin(s\phi)}{\sin \phi} \] Here \(\displaystyle \phi\) is the angle between the two unit quaternions: \(\displaystyle \cos(\phi) = q_0 \cdot q_1\).
Distance Metrics
There are multiple ways to calculate the distance between two quaternions. See Stackexchange.
- \(\displaystyle \theta = \cos^{-1}( \langle q_1, q_2 \rangle)\)
- Note that this loss can be misleading for rotations because of double cover:
- \(\displaystyle R(q_1) = R(-q_1)\) but \(\displaystyle \cos^{-1}(x) \neq \cos^{-1}(-x)\)
- Thus some people use the distance below.
- \(\displaystyle \theta = \cos^{-1}(2 \langle q_1, q_2 \rangle^2 - 1)\)
Spatial Alignment Problem
Also known as the RMSD problem.
Let \(\displaystyle \{y_k\}\) be a data array with N columns of D-dimensional points, known as the reference structure.
Let \(\displaystyle \{x_k\}\) be a data array of matched points, known as the test structure.
The goal is to rotation \(\displaystyle \{x_k\}\) by a rotation matrix \(\displaystyle R_{D} \in SO(D)\) to minimize the mean squared distance:
\[
\mathbf{S}_D = \sum_{k=1}^{N} \Vert R_D \cdot x_k - y_x \Vert^2
\]
The paper by Hanson focuses on 3D points.
Cross-term maximization
It can be shown that the least squares minimization problem can be converted to a cross-term maximization:
\(\displaystyle
\begin{aligned}
\mathbf{S}_D &= \sum_{k=1}^{N} \Vert R_D x_k - y_k \Vert^2\\
&= \sum_{k=1}^{N} (R_D x_k - y_x)^T(R_D x_k - y_k)\\
&= \sum_{k=1}^{N} ( x_k ^T R_D^T - y_k^T)(R_D x_k - y_k)\\
&= \sum_{k=1}^{N} ( x_k ^T R_D^T - y_k^T)(R_D x_k - y_k)\\
&= \sum_{k=1}^{N} - 2(R_D x_k) \cdot y_k\\
\end{aligned}
\)
So minimization of \(\displaystyle \mathbf{S}_D\) is equivalent to maximization of \(\displaystyle \sum_{k=1}^{N} (R_D x_k) \cdot y_k = tr R_d \cdot E\) where \(\displaystyle E_{ab} = \sum_{k=1}^{N} [x_k]_a [y_k]_b = [X \cdot Y^T]_{ab}\). This is also equivalent to solving \(\displaystyle \delta(q) = q M(E) q^t\).
Algebraic Solutions
There are expressions for solving the eigenvalues of \(\displaystyle M\).
Orientation-Frame Alignment Problem
An orientation frame are the columns of the \(\displaystyle SO(3)\) orientation matrix.
In this problem, we assume we have \(\displaystyle n\) paired orientation frames, \(\displaystyle \{p_k\}\) and \(\displaystyle \{r_k\}\) which we want to align using quaternion \(\displaystyle q\).
Geodesic arc-length distance
For two unit quaternions \(\displaystyle q_1\) and \(\displaystyle q_2\), the geodesic arc length is \(\displaystyle \alpha\) where \(\displaystyle q_1 \cdot q_2 = \cos \alpha\).
Because of sign ambiguity, \(\displaystyle R(q_2) = R(-q_2)\), we take the minimum value:
\[
d_{geodesic}(q_1, q_2) = \min(\alpha, \pi-\alpha) \text{ where }
0 \leq d_{geodesic}(q_1, q_2) \leq \frac{\pi}{2}
\]
The least-squares optimization is: \[ \begin{align*} \mathbf{S}_{geodesic} &= \sum_{k=1}^{N}(\operatorname{arccos} |(q*p_k) \cdot r_k|)^2\\ &= \sum_{k=1}^{N}(\operatorname{arccos} |q \cdot (r_k * \bar{p}_k)|)^2 \end{align*} \]
Chord distance
The geodesic least squares cannot be solved using linear algebra.
Instead you can use chord distance:
\(\displaystyle \begin{aligned} d_{\text{chord}}(q_1, q_2) &= \min(\Vert q_1 - q_2 \Vert, \Vert q_1 + q_2 \Vert)\\ & 0 \leq d_{\text{chord}}(q_1, q_2) \leq \sqrt{2} \end{aligned} \)
In this case, our least squares problem is: \[ \mathbf{S}_{\text{chord}} = \sum_{k=1}^{N}\left( \min( \Vert (q*p_k) - r_k \Vert, \Vert (q*p_k)+r_k \Vert ) \right)^2 \]
Combined Point + Frame Alignment Problem
Given a spatial profile matrix \(\displaystyle M(E)\) and orientation frame profile matrix \(\displaystyle U(S)\), we can normalize the scale to have unit-eigenvalue and then apply linear interpolation with dimensional constant \(\displaystyle \sigma\): \[ \Delta_{xf}(t, \sigma) = q \cdot \left[ (1-t) \frac{M(E)}{\epsilon_{x}} + t \sigma \frac{U(S)}{\epsilon_{f}} \right] \cdot q \]
Or alternatively, we can use slerp: \[ q(t) = \operatorname{slerp}(q_{x:opt}, q_{f:opt}, t) \]