新建網(wǎng)頁 1
從四元數(shù)轉(zhuǎn)換到矩陣
為了將角位移從四元數(shù)轉(zhuǎn)換到矩陣形式,可以利用旋轉(zhuǎn)矩陣,它能計算繞任意軸的旋轉(zhuǎn):

這個矩陣是用n和θ表示的,但四元數(shù)的分量是:
w = cos(θ/2)
x = nx sin(θ/2)
y = ny sin(θ/2)
z = nz sin(θ/2)
讓我們看看能否將矩陣變形以代入w、x、y、z,矩陣的9個元素都必須這樣做。幸運的是,這個矩陣的結(jié)構(gòu)非常好,一旦解出對角線上的一個元素,其他元素就能用同樣的方法求出。同樣,非對角線元素之間也是彼此類似的。
考慮矩陣對角線上的元素,我們將完整地解出m11,m22和m33解法與之類似:
m11 = nx2(1 - cosθ) + cosθ
我們將從上式的變形開始,變形方法看起來像是在繞圈子,但你馬上就能理解這樣做的目的:

現(xiàn)在需要消去cosθ項,而代之以包含cosθ/2或sinθ/2的項,因為四元數(shù)的元素都是用它們表示的,像以前那樣,設α=θ/2,先用α寫出cos的倍角公式,再代入θ:

上式是正確的,但它和其他的標準形式不同,即:
m11 = 1 - 2y2 - 2z2
實際上,還有其他的形式存在。最著名的一個形式是m11 = w2 + x2 - y2- z2,因為w2 + x2 + y2 + z2 = 1,所以這三種形式是等價的。現(xiàn)在回過頭來看看能不能直接導出其他標準形式,第一步,n是單位向量,nx2+ny2 + nz2 = 1,則1 - nx2 = ny2 + nz2。
m11 = 1 - (1 - nx2)(2sin2(θ/2))
= 1 - (ny2 +nz2)(2sin2(θ/2))
= 1 - 2ny2sin2(θ/2) - 2nz2sin2(θ/2)
= 1 - 2y2 - 2z2
元素m22和m33可以用同樣的方法求得。
讓我們來看看非對角線元素,它們比對角線元素簡單一些,以m12為例子:
m12 = nxny(1 - cosθ) + nzsinθ

其他非對角線元素可用同樣的方法導出。
最后,給出從四元素構(gòu)造的完整旋轉(zhuǎn)矩陣,如公式10.23所示:

從矩陣轉(zhuǎn)換到四元數(shù)
為了從旋轉(zhuǎn)矩陣中抽出相應的四元數(shù),可以直接利用公式 10.23,檢查對角線元素的和(也稱作矩陣的軌跡)得到:

通過使軌跡中三個元素中的兩個為負,可以用類似的方法求得其他三個元素:

不幸的是,這種方法并不總是能正確工作,因為平方根的結(jié)果總是正值。(更加準確地說,沒有選擇正根還是負根的依據(jù)。)但是,q和-q代表相同的方位,我們能任意選擇用非負根作為4個分量中的一個并仍能得到正確的四元數(shù),只是不能對四元數(shù)的所有4個數(shù)都用這種方法。
另一個技巧是檢查相對于對角線的對稱位置上元素的和與差:

那么應選用四種方法中的哪個呢?似乎最簡單的策略是總是先計算同一個分量,如w,然后再計算x、y、z。這伴隨著問題,如果w=0,除法就沒有定義;如果w非常小,將會出現(xiàn)數(shù)值不穩(wěn)定。Shoemake建議先判斷w、x、y、z中哪個最大(不用做平方根運算),根據(jù)上面的表,用矩陣對角線計算該元素,再用它計算其他三個。
下面的代碼用一種非常直接的方式實現(xiàn)了這個方法。
Listing 10.4: Converting a rotation matrix to a quaternion
// Input matrix:
float m11,m12,m13;
float m21,m22,m23;
float m31,m32,m33;
// Output quaternion
float w,x,y,z;
// Determine which of w, x, y, or z has the largest absolute value
float fourWSquaredMinus1 = m11 + m22 + m33;
float fourXSquaredMinus1 = m11 – m22 – m33;
float fourYSquaredMinus1 = m22 – m11 – m33;
float fourZSquaredMinus1 = m33 – m11 – m22;
int biggestIndex = 0;
float fourBiggestSquaredMinus1 = fourWSquaredMinus1;
if (fourXSquaredMinus1 > fourBiggestSquaredMinus1) {
fourBiggestSquaredMinus1 = fourXSquaredMinus1;
biggestIndex = 1;
}
if (fourYSquaredMinus1 > fourBiggestSquaredMinus1) {
fourBiggestSquaredMinus1 = fourYSquaredMinus1;
biggestIndex = 2;
}
if (fourZSquaredMinus1 > fourBiggestSquaredMinus1) {
fourBiggestSquaredMinus1 = fourZSquaredMinus1;
biggestIndex = 3;
}
// Perform square root and division
float biggestVal = sqrt(fourBiggestSquaredMinus1 + 1.0f) * 0.5f;
float mult = 0.25f / biggestVal;
// Apply table to compute quaternion values
switch (biggestIndex) {
case 0:
w = biggestVal;
x = (m23 – m32) * mult;
y = (m31 – m13) * mult;
z = (m12 – m21) * mult;
break;
case 1:
x = biggestVal;
w = (m23 – m32) * mult;
y = (m12 + m21) * mult;
z = (m31 + m13) * mult;
break;
case 2:
y = biggestVal;
w = (m31 – m13) * mult;
x = (m12 + m21) * mult;
z = (m23 + m32) * mult;
break;
case 3:
z = biggestVal;
w = (m12 – m21) * mult;
x = (m31 + m13) * mult;
y = (m23 + m32) * mult;
break;
}