How to construct a quaternion from two or three normal vectors that have direction and torque/roll of a solid body taken from a quaternion

Question

I am working with magnetic coils normal to each other (two or three) that a third party produces a quaternion that captures the direction and the roll around the direction. I must low pass filter the quaternions to smooth out noisy solid body orientation vibrations / movements. I am NOT asking how to filter.

After filtering (embedded code), I must reconstruct a new quaternion and send it along. However, the quaternion must behave the same as the third party.

This post uses angles to solve a problem and does not match the quaternion from the third party.

To obtain the direction of the solid body from the third party quaternion, I use $\mathbf {d'} =\mathbf {q} \mathbf { \hat{x} } \mathbf {q} ^{-1}$; I use obtain the normal for the attitude with $\mathbf {n'} =\mathbf {q} \mathbf { \hat{z} } \mathbf {q} ^{-1}$ where $\mathbf { \hat{x} }$ and $\mathbf { \hat{z} }$ are cartesian 3D basis vectors.

The quaternion has non-zero values for x, y, z, and w. There are plenty of examples of a rotation between two vectors. There appears scant information on preserving the rotation theta without using angles, as angles have gimble lock issues and or poles that are undefined. Quaternions do not have these issues. Vectors do not have these issues.

I am at a loss of how to reproduce the quaternion from these two filtered vectors: $\mathbf {d''}$ and $\mathbf {n''}$ where $\mathbf {n''}$ is normal to $\mathbf {d''}$ after filtering.

No matter how the quaternion $\mathbf {q'}$ is computed, it must have these properties: $\mathbf {d''} =\mathbf {q'} \mathbf { \hat{x} } \mathbf {q'} ^{-1}$ and $\mathbf {n''} =\mathbf {q'} \mathbf { \hat{z} } \mathbf {q'} ^{-1}$.

Given the above, how do I create $\mathbf {q'}$ from $\mathbf {d''}$ and $\mathbf {n''}$?

P.S.: I can get another norm from the quaternion and basis vector $\mathbf { \hat{y} }$, if this helps.

I am testing with this quaternion: $$ \begin{align*} q &= q_{w} + q_{x} \cdot i + q_{y} \cdot j + q_{z} \cdot k\\ q_{w} &= 3.666827976703644 \cdot 10^{-1}\\ q_{x} &= 8.160722255706787 \cdot 10^{-1}\\ q_{y} &= 3.313103318214417 \cdot 10^{-1}\\ q_{z} &= 2.996720373630524 \cdot 10^{-1}\\ \end{align*} $$

Answer 1

$ \renewcommand\vec\mathbf \newcommand\uvec[1]{\hat{\vec{#1}}} $We can do this by expanding the quaternions to the geometric algebra of 3D Euclidean space. This gives us a product on vectors, the geometric product, which represents reflections in the plane orthogonal to a vector $r$ as $v \mapsto -r^{-1}vr$. The product of two reflections $r_1r_2$ represents a rotation and is in fact a quaternion. The vectors themselves are not quaternions, but a vector $v$ can be converted to the corresponding quaternion $\vec v$ by $\vec v = -Iv$ where $I$ is a special element of the algebra called the right-handed unit pseudoscalar.

First reflect $\uvec x$ to $\vec d''$; this requires a vector $r_1$ (which we may take to be unit) such that $$ r_1^2 = 1,\quad d'' = -r_1\hat xr_1 $$ where e.g. $\hat x = I\uvec x$ since $I^{-1} = -I$. Having reflected everything with $r_1$, we then need to reflect the image $\hat z' = -r_1\hat zr_1$ of $\hat z$ to $n''$ while fixing $d''$; this requires using $r_2$ such that $$ r_2^2 = 1,\quad r_2\cdot d'' = 0,\quad n'' = -r_2\hat z'r_2. $$ We can find $\hat z'$ using standard vector algebra via $$ \hat z' = \hat z - 2(\hat z\cdot r_1)r_1 = \hat z - 2\frac{\hat z\cdot d''}{|d''-\hat x|^2}(d''-\hat x). $$ $r_1$ must be $(d''-\hat x)/|d''-\hat x|$ since both are unit vectors, similarly $r_2 = (n''-\hat z')/|n''-\hat z'|$. For consistency we then need $$ 0 = (n'' - \hat z')\cdot d'' = n''\cdot d''. $$ The last equality is because $d''$ is the image of $\hat x$ under $r_1$, so $\hat z'$ and $d''$ must be orthogonal.

Our rotation quaternion is $\vec q' = r_1r_2$; ignoring the normalization for a moment, $$ (d'' - \hat x)(n'' - \hat z') = (d'' - \hat x)\cdot(n'' - \hat z') + (d''-\hat x)\wedge(n''-\hat z'). $$ This wedge product is precisely the negative of the cross product of quaternion vectors, e.g. $v\wedge w = -\vec v\times\vec w$. To line up with your convention for rotating by a quaternion, we also need a final conjugation. All together, $$ \vec q' = \frac{(\vec d'' - \uvec x)\cdot(\vec n'' - \uvec z') + (\vec d'' - \uvec x)\times(\vec n'' - \uvec z')}{|\vec d''-\uvec x||\vec n'-\uvec z'|}. $$ It is probably better to ignore the denominator here and instead compute the numerator and normalize it.

---

Here is a rough Julia script demonstrating this algorithm. n is $n''$, d is $d''$, x is $\hat x$, z is $hat z$, z2 is $\hat z'$, and q2 is the quaternion determined by the algorithm. We start with a test quaternion q which we use to generate d and n from x and z, and then in the end we compare q2 with q using the norm of their difference. Note for any quaternion $q$ that $q$ and $-q$ generate the same rotation; we see here that we end up with q2 as the negative of q, which is why we compute norm2(q + q2) rather than norm2(q - q2).

using LinearAlgebra: cross, dot, norm2

qprod(q1, q2) = [q1[4]*q2[1:3] + q2[4]*q1[1:3] + cross(q1[1:3], q2[1:3]); q1[4]*q2[4] - dot(q1[1:3], q2[1:3])]
qconj(q) = [-q[1:3]; q[4]]
qinv(q) = qconj(q)/dot(q,q)
transform(q, v) = qprod(q, qprod(v, qconj(q)))

# x, y, z, w
q = [8.160722255706787e-01, 3.313103318214417e-01, 2.996720373630524e-01, 3.666827976703644e-01]
q = q/norm2(q)

x = [1, 0, 0]
z = [0, 0, 1]

d = transform(q, [x; 0])
@show d[4]
d = d[1:3]
d = d/norm2(d)

n = transform(q, [z; 0])
@show n[4]
n = n[1:3]
n = n/norm2(n)

@show dot(n, d)

z2 = z - 2dot(z, d)/dot(d-x, d-x)*(d-x)
z2 = z2/norm2(z2)
q2 = [cross(d-x, n-z2); dot(d-x, n-z2)]
q2 = q2/norm2(q2)

@show q
@show q2
@show norm2(q + q2)
```