Slerp

Hi all,

It was suggested that I do this a long time ago to compute the steps between two arbitrary rotations in an object animation. This code computes the rotation that is the specified fraction between its two argument rotations when measured along a geodesic curve in quaternion space. Practically this means that when going from arbitrary rot to another arbitrary rot you can easily find a rot that is 'between' them is a pleasing and intuitive sense.

Edit: I should point out that there are faster ways to do this, and I may even code one of them up in the near future. but this is working right now, in-world.

// -*- c++ -*-
// spherical linear interpolation library for LSL. cut and paste
// into your script! this will not compile on it's own (it has no
// default state definition)

float slerp_epsilon = 0.1;

float rot_dot(rotation l, rotation r) {
    return
        l.x * r.x
        + l.y * r.y
        + l.z * r.z
        + l.s * r.s
        ; }


rotation rot_neg(rotation r) {
    return
        <-r.x, -r.y, -r.z, -r.s>; }


float rot_mag(rotation r) {
    rotation q = <-r.x, -r.y, -r.z, r.s>;
    return
        llSqrt(r.x*q.x + r.y*q.y + r.z*q.z + r.s*q.s); }


rotation rot_norm(rotation r) {
    float mag = rot_mag(r);
    return
        <r.x / mag, r.y / mag, r.z / mag, r.s / mag>; }


// calculates the rotation that is *f* of the way between 
// *q1* and *q2*. *f* mudt be between 0 and 1
//
rotation rot_slerp(rotation q1, rotation q2, float f) {
    float r = 0;
    float s = 0;

    float lambda = rot_dot(q1, q2);
    if(lambda < 0) {
        q2 = rot_neg(q2);
        lambda = -lambda; }

    if(lambda < slerp_epsilon) {
        // linear interpolation
        r = 1 - f;
        s = f; }

      else {
        // spherical interpolation
        float alpha = llAcos(lambda);
        float gamma = 1 / llSin(alpha);
        r = llSin((1 - f) * alpha) * gamma;
        s = llSin(f * alpha) * gamma; }

    rotation q =
        <r * q1.x + s * q2.x,
         r * q1.y + s * q2.y,
         r * q1.z + s * q2.z,
         r * q1.s + s * q2.s>
        ;

    return rot_norm(q); }

An interesting bit of code for all three of us that know how to use it.

However, I am not seeing how this differs from simple quaternion division. Since a quaternion divided by another quaternion in Second Life equals the distance between the two, how does that differ from what you have here?

It looks to me like you're reinventing the wheel. Maybe I'm completely missing what it is you're doing here?



Edit: Aha. I see what's going on now. I'm using something similar, and simpler, with something I recently did for my game project.

Another method, albeit "cheap," is to just do this through Euler math. Simply convert the distance to a euler, divide, convert back.

1) the vector space you use to interpolate is non-euclidean, so while the distance is useful, it doesn't tell you how to interpolate in the same straghtforward way that you can with vectors.

2) this is interpolation, not distance. as a side effect the distance is calculated along the geodesic from quat to quat, not the straight-line distance.

(and no I don;t fully understand the algebraic geometry I'm spouting off here



I did something like that in an initial version and the LL rotation logic started to do funny things where my object would take the long way around between llSetRot() calls. Eggy Lippman first directed me towards SLERP as a 100% correct way to do this. I haven't had a bad surprise yet

- C

I think this code is equivalent.

rotation rotBetween( rotation a, rotation b, float f ) {
    float angleBetween = llAngleBetween(a, b);
    vector rotationalAxis = llRot2Axis(b/a);
    return a*llAxisAngle2Rot(rotationalAxis, angleBetween*f);
}

Whoops, there was a mistake in my post. Forgot to rotate the axis back, so it goes a little cock-eyed when A's natural axis of rotation wasn't parallel to the axis of rotation between a and b.

Fixed.

rotation rotBetween( rotation a, rotation b, float f ) {
    return a*llAxisAngle2Rot(llRot2Axis(b/a)*a, llAngleBetween(a, b)*f);
}

Where were you back in February? And where do I get a decent tutorial on quaternion algebra? Everything I've read so far leaves my head spinning...

- C

For what it's worth, I don't think it's really beneficial to understand the ins and outs of how quaternions work inside LSL. I think of rotations as just that - a data structure that describes a rotation. I don't much care that it has 4 elements or anything.

The important properties are:
- You can transform a vector with them
v*r = v'
- You can compose rotations together
r1 * r2 = r3
- Rotations are all invertible, and associative, but not commutative
1/r, (r1*r2)*r3 = r1*(r2*r3), r1*r2 != r2*r1
- You can use all these ll* functions to manipulate them.
http://secondlife.com/badgeo/wakka.php?wakka=Rotation

I haven't found a reason yet why anyone would need a deeper understanding than that inside LSL.

The only time where it becomes important to understand a quaternion is when you want to manipulate them outside LSL.

Heh heh... I disagree. But then, I'm an efficiency nazi with LSL and use these concepts all the time.

Two other excellent, "get your feet wet" guides on quaternion math:
http://en.wikipedia.org/wiki/Quaternion
http://www.cprogramming.com/tutorial/3d/quaternions.html


As well as an introduction to the linear algebra this is based on:
http://www.euclideanspace.com/maths/algebra/matrix/index.htm

Gah. Hit Quote instead of Edit. A+ to me.

From here:

I can't think of any time where manipulating quaternion elements directly is quicker or more compact than using the ll* API. Do you have any tricks to show us?

*necrobump*

Just in case anyone refers to this thread, I was playing around with this code today, and I realized that sometimes the "in-between" rotations took the "long way around", because llAngleBetween(a, b) likes returning positive values.

This one makes sure it always goes the short way

rotation rotBetween( rotation a, rotation b, float f ) {
    float angleBetween = llAngleBetween(a, b);
    
    if ( angleBetween > PI )
        angleBetween = angleBetween - TWO_PI;
    return a*llAxisAngle2Rot(llRot2Axis(b/a)*a, angleBetween*f);
}

[ necrobump due to:
/54/22/282969/1.html#post2164361
dated 20/09/2008 where this thread is referred to, and which provides the solution to the problem ]

Many thanks for this solution, and to Caoimhe Armitage for pointing it out

I'll have to test it, but I believe you should be able to get rid of the extra multiplication of the partial transformation axis by 'a' if you reverse the order of the final multiplication:

rotation interpolateRot(rotation a, rotation b, float t)
{
   rotation tAToB = b/a;  // Note that tAToB*a = (b/a)*a = b

   vector axis = llRot2Axis(tAToB);
   float angle = llRot2Angle(tAToB);
   if (angle > PI)
   {
      angle -= TWO_PI;
   }

   rotation tInterp = llAxisAngle2Rot(axis, t*angle); // = ZERO_ROT when t=0; tAToB when t=1

   return tInterp*a;  // = a when t=0; b when t=1
}

EDIT: Added back in the angle > PI fix above.

EDIT: Tested with a reasonable number of cases and it seems to work fine.

For those who want a short tight version...

rotation interpolateRot(rotation a, rotation b, float t)
{
   float angle = llRot2Angle(b /= a);
   if (angle > PI)
      angle -= TWO_PI;
   return llAxisAngle2Rot( llRot2Axis(b), t * angle) * a;
}

I've done a rundown on the math and it looks good.

I've been doing some digging through the documentation and the SL source.

Assuming the caveat on llRot2Angle is correct (which I trust that it is considering LL's math library does the same thing...) than you can eliminate the angle adjustment code. I have adjusted the reference implementations for both llRot2(Axis & Angle) to reflect this.

rotation slerp( rotation a, rotation b, float t ) {
   return llAxisAngle2Rot( llRot2Axis(b /= a), t * llRot2Angle(b)) * a;
}//Written collectively, Taken from http://forums.secondlife.com/showthread.php?p=536622

On a side note, because of the Axis sign changes used to ensuring that the angle is always positive and less than PI, llAngleBetween is a bad choice.

I would appreciate someone checking my reference implementations for llRot2Axis and llRot2Angle as at this time as i cannot get into SL. They need to be checked against the built in implementations.
https://wiki.secondlife.com/wiki/llRot2Axis
https://wiki.secondlife.com/wiki/llRot2Angle