一种另类的计算向量旋转公式（复杂慎用）

　　一般来说，我们解决向量旋转问题一般要么是用旋转矩阵，要么是用四元数。但很早以前我从网上找了一种比较另类的函数，当时也没有深究。最近又把这个函数拿出看看，仔细一琢磨，发现真的很另类。这里分享一下，就当是扩展一下思维。这个旋转方法据说叫Rodrigues旋转公式。

　　这个方法的公式是这样的，P'=P·cosθ + (A×P)sinθ +A(A·P)(1 - cosθ)。这种公式任谁第一眼看到都会摸不着头脑，我们首先来将公式换一个写法：

 图1

　　这个公式中|A|=1。为什么这样变换呢？因为P'正好是三个向量相加，后面会详细分析的，我们先将假设P'=u₁+u₂+u₃。接下来我们分析每一项具体是什么。如图2所示，P绕A轴旋转θ角，其中θ角如图位置。φ角是向量P和轴A之间的夹角。u₁是向量P在A轴的投影，u₂是向量(P'-u₁)在向量(P-u₁)的投影，u₃是从u₂末端到P'的向量，垂直于u₂。

 

图2

　　u₁的公式如下：

图3

　　其中|A|=1，u₂的公式如下：

图4

　　其中|P-u₁|=|P'-u₁|，u₃的公式如下：

图5

　　其中公式已经很清楚的阐述了推理过程，下面贴出实现的函数代码。

 

void rotateWithAxis(CVector3& rotate, CVector3& axis, GLdouble angle, CVector3& direct)
{
	/*
		P绕A轴旋转角度θ，得到新的向量P'，则：
		P'=P * cosθ + (A×P)sinθ +A(A·P)(1 - cosθ)
		其中A为单位向量，旋转角度θ为"逆时针方向"旋转的角度。

		rotate 表示的是向量P
		axis     表示的是轴A
		angle   表示的是角度θ
		direct  表示的是向量P'

		A×P=（ay * pz- az * py, ax * pz- az * px , ax * py- ay * px）
		A·P = ax * px + ay * py + az * pz
	*/
#define PI 3.1415926535898
#define Cos cos(anglePI)
#define Sin sin(anglePI)
	GLdouble anglePI = angle/180 * PI;
	CVector3 cross;
	GLdouble s;
	axis.Normal();
	cross = CVector3::Cross(rotate, axis);
	s = rotate[0]*axis[0] + rotate[1]*axis[1] + rotate[2]*axis[2];
	
	direct.SetVector3(rotate[0]*Cos + cross[0]*Sin + axis[0]*s*(1-Cos),
		rotate[1]*Cos + cross[1]*Sin + axis[1]*s*(1-Cos), rotate[2]*Cos + cross[2]*Sin + axis[2]*s*(1-Cos));
}