How would I apply an Exponential Moving Average to Quaternions? I'm trying to filter positional and rotational data using an Exponential Moving Average (EMA) filter. This has worked without issues for positional data (3D vectors) but I can't figure it out for rotational data (Quaternions). I first tried the naive way, which is applying the EMA to each separate component, but this obviously doesn't work very well. How would I do it?
 A: Equations 12 and 13 in Markley et al.'s paper "Quaternion Averaging" (2007) define the weighted average of $n$ quaterions $\mathbf q_i$ with scalar weights $w_i$:
$$\bar{\mathbf q}=\arg\,\max_{\mathbf q\in\mathbb S^3}\mathbf q^T\mathbf M\mathbf q$$
where
$$\mathbf M=\sum_{i=1}^nw_i\mathbf q_i\mathbf q_i^T.$$
That is, interpret the quaternions as vectors in $\mathbb R^4$, form the matrix $\mathbf M$, and take its largest eigenvector.
Here, by the "average quaternion" we mean the one that minimizes the weighted sum of squared differences of the corresponding rotation matrices, measured in the Frobenius norm (equation 2 of the above paper).
A: One of the standard techniques in using quaternions for graphics applications is the notion of 'slerping', Spherical Linear Interpolation, for interpolating between two quaternions.  The technique is an extension of the standard form of linear interpolation to work along spherical arcs ('moving' at constant angular velocity): essentially, rather than computing $(1-\alpha)v_0+\alpha v_1$, you instead compute $\dfrac{\sin((1-\alpha)\theta)}{\sin\theta}\mathbf{q}_0+\dfrac{\sin(\alpha\theta)}{\sin\theta}\mathbf{q}_1$, where $\cos\theta=\mathbf{q}_0\cdot\mathbf{q}_1$ (i.e., $\theta$ is the angle between the two vectors). It's a nice exercise to prove that this produces a unit quaternion, and that it reduces to linear interpolation in the limit of small $\theta$.
Since every step of the EMA filter is an interpolation between the previous smoothed value and the new value, you can use slerping to perform those interpolations on your quaternion data.  The operation is somewhat expensive for real-time calculation (it involves an inverse trigonometric function evaluation and multiple trig function evaluations), but it's so fundamental that many articles have been written on making faster approximate versions, and searching for 'optimized slerp' or 'fast slerp' should yield plenty of good results for you.
