这是国内找不到的超好文章。（为什么大陆的4元数文章很垃圾呢？）（翻译中。。。奉献给大家~~）
70秒即懂，能使用，用四元数，4元数，阔特尼恩，Quaternion旋转(C) 中田 亨 （独立行政法人 产业技术综合研究所 数字人类研究中心 研究员 博士（工学）） 2003年11月25日
★这个页面的对象读者想把三次元的旋转，用CG等定量地处理的人使用欧拉角(Euler Angles)的话，不懂得其道理的人 卡尔丹角和欧拉角(Cardan Angles)不能区别的人 对吉恩瓦尔洛克很困惑的人但是，对数学之类麻烦的事情很讨厌的人想要实例程序的人 没有时间的人
★旋转篇： 我将说明使用了四元数（si yuan shu, quaternion）的旋转的操作步骤（１）四元数的虚部，实部和写法所谓四元数，就是把4个实数组合起来的东西。4个元素中，一个是实部，其余3个是虚部。比如，叫做Q的四元数，实部t而虚部是x,y,z构成，则像下面这样写。Q = (t; x, y, z) 又，使用向量 V=(x,y,z)，Q = (t; V) 也可以这么写。
正规地用虚数单位i,j,k的写法的话，Q = t + xi + yj + zk 也这样写，不过，我不大使用
（２）四元数之间的乘法虚数单位之间的乘法 ii = -1, ij = -ji = k (其他的组合也是循环地以下同文) 有这么一种规则。（我总觉得，这就像是向量积（外积），对吧） 用这个规则一点点地计算很麻烦，所以请用像下面这样的公式计算。
A = (a; U) B = (b; V) AB = (ab - U·V; aV + bU + U×V)不过，“U·V”是内积，「U×V」是外积的意思。注意：一般AB<>BA所以乘法的左右要注意！
（3）3次元的坐标的四元数表示如要将某坐标(x,y,z)用四元数表示，P = (0; x, y, z) 则要这么写。 另外，即使实部是零以外的值，下文的结果也一样。用零的话省事所以我推荐。
（４）旋转的四元数表示以原点为旋转中心，旋转的轴是(α, β, γ)（但 α^2 + β^2 + γ^2 = 1）， （右手系的坐标定义的话，望向向量(α, β, γ)的前进方向反时针地） 转θ角的旋转，用四元数表示就是，Q = (cos(θ/2); α sin(θ/2), β sin(θ/2), γ sin(θ/2)) R = (cos(θ/2); -α sin(θ/2), -β sin(θ/2), -γ sin(θ/2)) (另外R 叫 Q 的共轭四元数。）
那么，如要实行旋转，则 R P Q = (0; 答案)
请像这样三明治式地计算。这个值的虚部就是旋转之后的点的坐标值。 （另外，实部应该为零。请验算看看）
/// Quaternion.cpp /// (C) Toru Nakata, toru-nakata@aist.go.jp /// 2004 Dec 29 #include <math.h> #include <iostream.h> /// Define Data type typedef struct { double t; // real-component double x; // x-component double y; // y-component double z; // z-component } quaternion; Kakezan quaternion Kakezan(quaternion left, quaternion right) { quaternion ans; double d1, d2, d3, d4;
d1 = left.t * right.t; d2 = -left.x * right.x; d3 = -left.y * right.y; d4 = -left.z * right.z; ans.t = d1+ d2+ d3+ d4;
d1 = left.t * right.x; d2 = right.t * left.x; d3 = left.y * right.z; d4 = -left.z * right.y; ans.x = d1+ d2+ d3+ d4;
d1 = left.t * right.y; d2 = right.t * left.y; d3 = left.z * right.x; d4 = -left.x * right.z; ans.y = d1+ d2+ d3+ d4;
d1 = left.t * right.z; d2 = right.t * left.z; d3 = left.x * right.y; d4 = -left.y * right.x; ans.z = d1+ d2+ d3+ d4; return ans; } Make Rotational quaternion quaternion MakeRotationalQuaternion(double radian, double AxisX, double AxisY, double AxisZ) { quaternion ans; double norm; double ccc, sss;
ans.t = ans.x = ans.y = ans.z = 0.0;
norm = AxisX * AxisX + AxisY * AxisY + AxisZ * AxisZ; if(norm <= 0.0) return ans;
norm = 1.0 / sqrt(norm); AxisX *= norm; AxisY *= norm; AxisZ *= norm;
ccc = cos(0.5 * radian); sss = sin(0.5 * radian);
ans.t = ccc; ans.x = sss * AxisX; ans.y = sss * AxisY; ans.z = sss * AxisZ;
return ans; } Put XYZ into quaternion quaternion PutXYZToQuaternion(double PosX, double PosY, double PosZ) { quaternion ans;
ans.t = 0.0; ans.x = PosX; ans.y = PosY; ans.z = PosZ;
return ans; } / main int main() { double px, py, pz; double ax, ay, az, th; quaternion ppp, qqq, rrr;
cout << "Point Position (x, y, z) " << endl; cout << " x = "; cin >> px; cout << " y = "; cin >> py; cout << " z = "; cin >> pz; ppp = PutXYZToQuaternion(px, py, pz);
while(1) { cout << "/nRotation Degree ? (Enter 0 to Quit) " << endl; cout << " angle = "; cin >> th; if(th == 0.0) break;
cout << "Rotation Axis Direction ? (x, y, z) " << endl; cout << " x = "; cin >> ax; cout << " y = "; cin >> ay; cout << " z = "; cin >> az;
th *= 3.1415926535897932384626433832795 / 180.0; /// Degree -> radian;
qqq = MakeRotationalQuaternion(th, ax, ay, az); rrr = MakeRotationalQuaternion(-th, ax, ay, az);
ppp = Kakezan(rrr, ppp); ppp = Kakezan(ppp, qqq);
cout << "/nAnser X = " << ppp.x << "/n Y = " << ppp.y << "/n Z = " << ppp.z << endl;
}
return 0; }
