这是国内找不到的超好文章。（为什么大陆的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; 
}  

/// 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; 
}  

/// 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; 
}  

恩，支持