今天我们和大家分享的是四元数的插值。这里的插值指的是球面线性插值。例如,我们要模拟一下地球绕着太阳,从P1到P2。这中间的每一个位置都要用球面线性插值来做。
首先我们聊一聊线性插值:
x=x1-x2,t是插值系数,则lerp(x1,x2,t)=x1+t*x表示x1到x2的插值。
四元数的插值
q=(p-1)*p1,插值系数为t,则p到p1的插值为:slerp(p,p1,t)=p*((p-1)p1)t,表示p的逆乘以p1,他们乘积的t次方,乘以p
但是上面的公式在编程的时候用的挺麻烦的,所以我们使用下面的公式:
旋转插值,想开头所说的,这里说的插值是球面插值,是在3D空间中旋转的,因此我们可以将它等价于旋转插值。
例如:向量V1,V0
W=v1-v0
vt=v0+tw
通过上图可知,VT=K0V0+K1V1.在这里,V1,V0,VT都是单位向量,V1和K1V1平行。
因此,我们可以求得sinw=sintw/k1,则k1=sintw/sinw,同样可得 k0=sin(1-t)w/sinw.
所以代入VT,就可以求出VT了。
四元数也是一个道理,注意slerp这个函数返回的是一个四元数。
以上就是四元数插值的公式,在编程实现的时候我们要注意
使用点成来纠结夹角W,还有就是放夹角W很小的时候,sinw很小,但是cosw趋于1,因此,就变成了线性插值。
在这里,我们贴出cocos2dx 3.6中四元数的实现代码,大家学习学习。
/** * Interpolates between two quaternions using spherical linear interpolation. * * Spherical linear interpolation provides smooth transitions between different * orientations and is often useful for animating models or cameras in 3D. * * Note: For accurate interpolation,the input quaternions must be at (or close to) unit length. * This method does not automatically normalize the input quaternions,so it is up to the * caller to ensure they call normalize beforehand,if necessary. * * @param q1x The x component of the first quaternion. * @param q1y The y component of the first quaternion. * @param q1z The z component of the first quaternion. * @param q1w The w component of the first quaternion. * @param q2x The x component of the second quaternion. * @param q2y The y component of the second quaternion. * @param q2z The z component of the second quaternion. * @param q2w The w component of the second quaternion. * @param t The interpolation coefficient. * @param dstx A pointer to store the x component of the slerp in. * @param dsty A pointer to store the y component of the slerp in. * @param dstz A pointer to store the z component of the slerp in. * @param dstw A pointer to store the w component of the slerp in. */ static void slerp(float q1x,float q1y,float q1z,float q1w,float q2x,float q2y,float q2z,float q2w,float t,float* dstx,float* dsty,float* dstz,float* dstw);
void Quaternion::slerp(float q1x,float* dstw) { // Fast slerp implementation by kwhatmough: // It contains no division operations,no trig,no inverse trig // and no sqrt. Not only does this code tolerate small constraint // errors in the input quaternions,it actually corrects for them. GP_ASSERT(dstx && dsty && dstz && dstw); GP_ASSERT(!(t < 0.0f || t > 1.0f)); if (t == 0.0f) { *dstx = q1x; *dsty = q1y; *dstz = q1z; *dstw = q1w; return; } else if (t == 1.0f) { *dstx = q2x; *dsty = q2y; *dstz = q2z; *dstw = q2w; return; } if (q1x == q2x && q1y == q2y && q1z == q2z && q1w == q2w) { *dstx = q1x; *dsty = q1y; *dstz = q1z; *dstw = q1w; return; } float halfY,alpha,beta; float u,f1,f2a,f2b; float ratio1,ratio2; float halfSecHalfTheta,versHalfTheta; float sqNotU,sqU; float cosTheta = q1w * q2w + q1x * q2x + q1y * q2y + q1z * q2z; // As usual in all slerp implementations,we fold theta. alpha = cosTheta >= 0 ? 1.0f : -1.0f; halfY = 1.0f + alpha * cosTheta; // Here we bisect the interval,so we need to fold t as well. f2b = t - 0.5f; u = f2b >= 0 ? f2b : -f2b; f2a = u - f2b; f2b += u; u += u; f1 = 1.0f - u; // One iteration of Newton to get 1-cos(theta / 2) to good accuracy. halfSecHalfTheta = 1.09f - (0.476537f - 0.0903321f * halfY) * halfY; halfSecHalfTheta *= 1.5f - halfY * halfSecHalfTheta * halfSecHalfTheta; versHalfTheta = 1.0f - halfY * halfSecHalfTheta; // Evaluate series expansions of the coefficients. sqNotU = f1 * f1; ratio2 = 0.0000440917108f * versHalfTheta; ratio1 = -0.00158730159f + (sqNotU - 16.0f) * ratio2; ratio1 = 0.0333333333f + ratio1 * (sqNotU - 9.0f) * versHalfTheta; ratio1 = -0.333333333f + ratio1 * (sqNotU - 4.0f) * versHalfTheta; ratio1 = 1.0f + ratio1 * (sqNotU - 1.0f) * versHalfTheta; sqU = u * u; ratio2 = -0.00158730159f + (sqU - 16.0f) * ratio2; ratio2 = 0.0333333333f + ratio2 * (sqU - 9.0f) * versHalfTheta; ratio2 = -0.333333333f + ratio2 * (sqU - 4.0f) * versHalfTheta; ratio2 = 1.0f + ratio2 * (sqU - 1.0f) * versHalfTheta; // Perform the bisection and resolve the folding done earlier. f1 *= ratio1 * halfSecHalfTheta; f2a *= ratio2; f2b *= ratio2; alpha *= f1 + f2a; beta = f1 + f2b; // Apply final coefficients to a and b as usual. float w = alpha * q1w + beta * q2w; float x = alpha * q1x + beta * q2x; float y = alpha * q1y + beta * q2y; float z = alpha * q1z + beta * q2z; // This final adjustment to the quaternion's length corrects for // any small constraint error in the inputs q1 and q2 But as you // can see,it comes at the cost of 9 additional multiplication // operations. If this error-correcting feature is not required,// the following code may be removed. f1 = 1.5f - 0.5f * (w * w + x * x + y * y + z * z); *dstw = w * f1; *dstx = x * f1; *dsty = y * f1; *dstz = z * f1; }