Component of a quaternion rotation around an axis

minorlogic's answer pointing out the swing-twist decomposition is the best answer so far, but it is missing an important step. The issue (using minorlogic's pseudocode symbols) is that if the dot product of ra and p is negative, then you need to negate all four components of twist, so that the resulting rotation axis points in the same direction as direction. Otherwise if you are trying to measure the angle of rotation (ignoring the axis), you'll confusingly get a mix of correct rotations and the reverse of the correct rotations, depending on whether or not the rotation axis direction happened to flip when calling projection(ra, direction). Note that projection(ra, direction) computes the dot product, so you should reuse it and not compute it twice.

Here's my own version of the swing-twist projection (using different variable names in some cases instead of minorlogic's variable names), with the dot product correction in place. Code is for the JOML JDK library, e.g. v.mul(a, new Vector3d()) computes a * v, and stores it in a new vector, which is then returned.

/**
 * Use the swing-twist decomposition to get the component of a rotation
 * around the given axis.
 *
 * N.B. assumes direction is normalized (to save work in calculating projection).
 * 
 * @param rotation  The rotation.
 * @param direction The axis.
 * @return The component of rotation about the axis.
 */
private static Quaterniond getRotationComponentAboutAxis(
            Quaterniond rotation, Vector3d direction) {
    Vector3d rotationAxis = new Vector3d(rotation.x, rotation.y, rotation.z);
    double dotProd = direction.dot(rotationAxis);
    // Shortcut calculation of `projection` requires `direction` to be normalized
    Vector3d projection = direction.mul(dotProd, new Vector3d());
    Quaterniond twist = new Quaterniond(
            projection.x, projection.y, projection.z, rotation.w).normalize();
    if (dotProd < 0.0) {
        // Ensure `twist` points towards `direction`
        twist.x = -twist.x;
        twist.y = -twist.y;
        twist.z = -twist.z;
        twist.w = -twist.w;
        // Rotation angle `twist.angle()` is now reliable
    }
    return twist;
}

I tried to implement sebf's answer, it seems good, except that the choice of the choice of vector in step 1:

1. Take the axis you want to find the rotation around, and find an orthogonal vector to it.

is not sufficient for repeatable results. I have developed this on paper, and I suggest the following course of action for the choice of the vector orthogonal to the "axis you want to find the rotation around", i.e. axis of observation. There is a plane orthogonal to the axis of observation. You have to project the axis of rotation of your quaternion onto this plane. Using this resulting vector as the vector orthogonal to the axis of observation will give good results.

Thanks to sebf for setting me down the right course.

Code for Unity3d

// We have some given data
Quaternion rotation = ...;
Vector3 directionAxis = ...;

// Transform quaternion to angle-axis form
rotation.ToAngleAxis(out float angle, out Vector3 rotationAxis);

// Projection magnitude is what we found - a component of a quaternion rotation around an axis to some direction axis
float proj = Vector3.Project(rotationAxis.normalized, directionAxis.normalized).magnitude;