#32 Bug in cgtypes.slerp() or cgtypes.quat.normalize()


There seems to be a bug in either cgtypes.slerp() or cgtypes.quat.normalize()

The documentation for slerp() states that "q0 and q1 must be unit quaternions". If they're not, the result is always (-1.#IND, -1.#IND, -1.#IND, -1.#IND)

So that's fine. I normalize my quats before passing them to slerp(). However, normalize() apparently doesn't always return a perfectly unit quaternion. For instance, I sometimes see a result like this: (1, -6.93889e-018, -8.67362e-019, -5e-007). And passing this to slerp() then causes the issue stated above.

I've tried to boil my issue down to a minimal example, by using some print statements to get the "offensive" quat values before and after normalization but there seems to be some inaccuracy in the way floats are displayed(?), and the bug is a little elusive.

At any rate, this example should make my confusion apparent:

weight = 0.14285714285714285
quatA = cgtypes.quat(1, -6.93889e-018, -8.67362e-019, -5e-007)
quatB = cgtypes.quat(1, 2.08167e-017, -3.98986e-017, -5e-007)
quatANorm, quatBNorm = quatA.normalize(), quatB.normalize()
quatC = cgtypes.slerp(weight, quatA, quatB)
quatCNorm = cgtypes.slerp(weight, quatANorm, quatBNorm)
print weight, quatA, quatB, quatC
print weight, quatANorm, quatBNorm, quatCNorm

0.142857142857 (1, -6.93889e-018, -8.67362e-019, -5e-007) (1, 2.08167e-017, -3.98986e-017, -5e-007) (-1.#IND, -1.#IND, -1.#IND, -1.#IND)
0.142857142857 (1, -6.93889e-018, -8.67362e-019, -5e-007) (1, 2.08167e-017, -3.98986e-017, -5e-007) (1, -6.93889e-018, -8.67362e-019, -5e-007)

As you can see, the values before / after normalization are the same (and not perfectly normalized!). Yet slerp() seems to know there is a difference, and returns a valid quat for the normalized ones.

Any ideas? This might just be due to my relative inexperience with Python...



  • Comment has been marked as spam. 

    You can see all pending comments posted by this user  here


    Anonymous - 2011-07-26

    After some more investigation, most other libraries I've looked at handle this issue (of very similar quats, or non-unit quats) in the slerp() method. And they just accept that normalize() may not return truly unit quats.

    See http://cgkit.svn.sourceforge.net/viewvc/cgkit/cgkit/trunk/supportlib/include/quat.h?revision=341&content-type=text%2Fplain (Line 187)

    acos(ca) results in a complex number for values of ca > 1 or ca < -1. I've rewritten my slerp function as follows to handle this better:

    #This is a Python version (but the same would work for C)
    def slerp(t,q0,q1,shortest=True):
    neg_q1 = False #Does q1 have to be negated (so that the shortest path is taken)?

    ca = q0.dot(q1)
    if (shortest and ca < 0):
    ca = -ca
    neg_q1 = True

    if ( (1.0 - ca) > 1e-12): #epsilion?
    #standard case
    o = math.acos(ca)
    so = math.sin(o)
    a = math.sin(o*(1.0-t)) / so
    b = math.sin(o*t) / so
    #quats are very close, do linear interpolation
    a = 1.0 - t
    b = t

    #do the interpolation
    if (neg_q1):
    return q0*a - q1*b
    return q0*a + q1*b

  • Matthias Baas

    Matthias Baas - 2011-07-26

    quatA and quatB in your example are not normalized whereas quatANorm and quatBNorm are which is why slerp() fails on the former but succeeds on the latter.
    When you print the quaternions you only see the actual value up to a certain precision and as the normalized quats are fairly close to the input quats, you don't see the difference in your prints.

    You can print the lengths or the individual components to see the difference:

    >>> abs(quatA)
    >>> abs(quatANorm)

    >>> quatA.x
    >>> quatANorm.x

    As I see in your second comment you got down to the bottom of it yourself and it's the acos() that is not defined for values >1 or <-1. Even though the above example doesn't really demonstrate the problem (as the input quats are clearly not normalized), it is not impossible that the length of the quaternion is still a slightly above 1 after normalization (because of numerical inaccuracies). So I'll change the code accordingly and make sure acos() won't be called with invalid values.
    Thanks for reporting this!

  • Matthias Baas

    Matthias Baas - 2011-08-16

    Fixed and committed.

  • Matthias Baas

    Matthias Baas - 2011-08-16
    • status: open --> closed-fixed

Log in to post a comment.