基矩阵 - 王朝网络

一、空间坐标系的基和基矩阵

在3-D空间中，我们用空间坐标系来规范物体的位置，空间坐标系由3个相互垂直的坐标轴组成，我们就把它们作为我们观察3-D空间的基础，空间中物体的位置可以通过它们来衡量。当我们把这3个坐标轴上单位长度的向量记为3个相互正交的单位向量i,j,k，空间中每一个点的位置都可以被这3个向量线性表出，如P<1,-2,3>这个点可以表为i-2j+3k。

我们把这3个正交的单位向量称为空间坐标系的基，它们单位长度为1且正交，所以可以成为标准正交基。三个向量叫做基向量。现在我们用矩阵形式写出基向量和基。

i = | 1 0 0 |

j = | 0 1 0 |

k = | 0 0 1 |

B = | i | | 1 0 0 |

| j | | 0 1 0 |

| k | | 0 0 1 |

这样的矩阵我们叫它基矩阵。有了基矩阵，我们就可以把空间坐标系中的一个向量写成坐标乘上基矩阵的形式，比如上面的向量P可以写成：

P = C x B=> | 1 0 0 | | 1 -2 3 | = | 1 -2 3 | x | 0 1 0 | | 0 0 1 |

这样的话，空间坐标系下的同一个向量在不同的基下的坐标是不同的。

二、局部坐标系和局部坐标

和空间坐标系（也可以叫做全局坐标系或者世界坐标系）并存的称为局部坐标系（也叫坐标架——coordinate frame），它有自己的基，这些基向量把空间坐标系作为参考系。比如

| x'| | -1 0 0 | B' = | y'| = | 0 1 0 | | z'| | 0 0 -1 | | x''| | 2^½ /2 0 2^½ /2 |

B'' = | y''| = | 0 -1 0 |

| z''| | -(2^½) /2 0 2^½ /2 |就是两个局部坐标系的基，如图:

现在我们可以把上面那个空间坐标中的向量P|1 -2 3|（以后都用矩阵表示）表示在不同的基下，我把它写成一个大长串的式子： | x' | | x''| P = | Px' Py' Pz' | x | y' | = | Px'' Py'' Pz'' | x | y''|

| z' | | z''|

这里| Px' Py' Pz'|是P在B'下的坐标，| Px'' Py'' Pz''|是P在B''下的坐标，我把它写的具体点吧：

| -1 0 0 | | 2^½ /2 0 2^½ /2| | 1 -2 3 | = | -1 -2 -3 | x | 0 1 0 | = | 2*2^½ -2 2^½ | x | 0 -1 0 |

| 0 0 -1 | | -(2^½) /2 0 2^½ /2|

这就是说，在空间坐标系下面的向量| 1 -2 3 |在基B'下的坐标为|-1 -2 -3|，在B''下的坐标为| 2*2^½ -2 2^½ |。当然空间坐标系也有自己的基B|i j k|^T（因为是列向量，所以写成行向量的转置），但我们现在是拿它当作一个参考系。

在研究了局部坐标系之后，我现在要分析两个应用它们的例子，先来看

三、空间坐标系中一个点围绕任一轴的旋转

上一篇讨论3-D空间旋转的时候说到有一个高档的方法做3-D空间任意轴旋转，现在我们的知识储备已经足够理解这个方法了(Quake引擎使用的就是这个方法)。

如上所示，空间坐标系中的一个局部坐标系xyz中有一个向量a(2,5,3)和一个点p(8,4,2)现在我要让p点围绕a向量旋转60度，得到p’点，该如何做呢？从目前掌握的旋转知识来看，我们有两个理论基础：

1）在一个坐标系中的一个点，如果要它围绕该坐标系中一个坐标轴旋转，就给它的坐标值乘相应的旋转矩阵，如

[cosA -sinA 0 ][sinA cosA 0 ][0 0 1 ]

等等。

2）我们已经学习了局部坐标系的理论了，知道空间中一个点在不同的坐标系中的坐标不同。利用这一点，我们可以很方便的让一个点或者向量在不同的坐标系之间转换。

我们联系这两个理论根据，得出我们的思路：

1构造另一个局部坐标系abc，使得a成为该坐标系的一个坐标轴。

2 把p的坐标变换到abc中，得到p’，用旋转公式让p’围绕已经成为坐标轴的a旋转，得到p’’。

3把p’’再变换回坐标系xyz，得到p’’’，则p’’’就是p围绕a旋转后的点。

下面我们逐步说明。

首先我们构造abc，我们有无数种方法构造，因为只要保证b、c之间以及他们和a之间都正交就可以了，但我们只要一个。根据上图，我们首先产生一个和a正交的b。这可以通过向量的叉乘来完成：我们取另一个向量v（显然，这个向量是不能和a共线的任何非零向量），让它和a决定一个平面x，然后让v叉乘a得到一个垂直于x的向量b，因为b垂直于x，而a在平面x上，因此b一定垂直于a，然后用a叉乘b得到c，最后单位化a、b、c，这样就得到了局部坐标系abc。

然后我们把p点变换到abc坐标系中，得到p’，即p’就是p在abc中的坐标：

|a b c| * p’= |x y z| * p

p’ = |a b c|^-1 * |x y z| * p

|ax bx cx| |1 0 0| |px|

p’ = |ay by cy| ^-1 * |0 1 0| * |py|

|az bz cz| |0 0 1| |pz|

注意这里|a b c|^-1即矩阵|a b c|的逆矩阵，因为a、b、c是三个正交向量，并且是单位向量，因此|a b c|是一个正交矩阵，正交矩阵的转置和逆相等，这是它的一个特性，因此上面的公式就可以写成:

|ax ay az| |1 0 0| |px|

p’ = |bx by bz| * |0 1 0| * |py|

|cx cy cz| |0 0 1| |pz|

这个时候p’就是p在abc坐标系下的坐标了。此时a已经是一个坐标轴了，我们可以用旋转矩阵来做。

p’’ = RotMatrix * p’

[1 0 0] |p’x|p’’ = [0 cos60 -sin60] * |p’y| [0 sin60 cos60] |p’z|

最后，我们把p’’再次变换回xyz坐标系，得到最终的p’’’

|a b c| * p’’ = |x y z| * p’’’

p’’’ = |x y z|^-1 * |a b c| * p’’

p’’’ = |a b c| * p’’

最后

p’’’ = |a b c| * RotMatrix * |a b c|^T * p = M * p

这样就得到了xyz坐标系中点p围绕a旋转60度后的点。

最后，我用Quake3引擎的相应函数(来自idSoftware ——quake3-1[1].32b-source——mathlib.c)来完成对这个算法的说明：

===============

RotatePointAroundVector

dst是一个float[3]，也就是p’’’

dir相当于a，point就是p，degrees是旋转度数

===============

void RotatePointAroundVector( vec3_t dst, const vec3_t dir, const vec3_t point,

float degrees ) {

float m[3][3];

float im[3][3];

float zrot[3][3];

float tmpmat[3][3];

float rot[3][3];

int i;

vec3_t vr, vup, vf;

float rad;

vf[0] = dir[0];

vf[1] = dir[1];

vf[2] = dir[2];

// 首先通过dir得到一个和它垂直的vr

// PerpendicularVector()函数用于构造和dir垂直的向量

// 也就是我们上面的第1步

PerpendicularVector( vr, dir );

// 通过cross multiply得到vup

// 现在已经构造出坐标轴向量vr, vup, vf

CrossProduct( vr, vf, vup );

// 把这三个单位向量放入矩阵中

m[0][0] = vr[0];

m[1][0] = vr[1];

m[2][0] = vr[2];

m[0][1] = vup[0];

m[1][1] = vup[1];

m[2][1] = vup[2];

m[0][2] = vf[0];

m[1][2] = vf[1];

m[2][2] = vf[2];

// 产生转置矩阵im

memcpy( im, m, sizeof( im ) );

im[0][1] = m[1][0];

im[0][2] = m[2][0];

im[1][0] = m[0][1];

im[1][2] = m[2][1];

im[2][0] = m[0][2];

im[2][1] = m[1][2];

// 构造旋转矩阵zrot

memset( zrot, 0, sizeof( zrot ) );

zrot[0][0] = zrot[1][1] = zrot[2][2] = 1.0F;

rad = DEG2RAD( degrees );

zrot[0][0] = cos( rad );

zrot[0][1] = sin( rad );

zrot[1][0] = -sin( rad );

zrot[1][1] = cos( rad );

// 开始构造变换矩阵M

// tmpmat = m * zrot

MatrixMultiply( m, zrot, tmpmat );

// rot = m * zrot * im

MatrixMultiply( tmpmat, im, rot );

// 则 rot = m * zrot * im 和我们上面推出的

// M = |a b c| * RotMatrix * |a b c|^T 一致

// 变换point这个点

// p’’’ = M * p

for ( i = 0; i < 3; i++ ) {

dst[i] = rot[i][0] * point[0] + rot[i][1] * point[1] + rot[i][2] * point[2];

}

四、世界空间到相机空间的变换

空间坐标系XYZ，相机坐标系UVN。这时候相机空间的基（以下简称相机）在空间坐标系中围绕各个坐标轴旋转了一定角度<a,b,c>，然后移动了<x,y,z>。对于模型我们可以看作相对于相机的逆运动，即模型旋转了一定角度<-a,-b,-c>，然后移动了<-x,-y,-z>，可以把相机和物体的运动看成两个互逆的变换。这样，可以通过对相机的变换矩阵求逆来得到模型的变换矩阵。下面来具体看一下，如何得到相机变换矩阵，并且求得它的逆矩阵。

首先声明一下，对于一个模型的变换，我们可以给模型矩阵左乘变换矩阵：

M x P = P'

| E F G H | | y | | Ex + Fy + Gz + H | x = | I J K L | | z | | Ix + Jy + Kz + L |

也可以右乘变换矩阵：

P^T x M^T = P'^T

| A E I M |

| B F J N || x y z 1| x = |Ax+By+Cz+D Ex+Fy+Gz+H Ix+Jy+Kz+L Mx+Ny+Oz+P| | C G K O |

| D H L P |

可以看出两种变换方式是一个转置关系，结果只是形式上的不同，但这里我们使用后者，即右乘变换矩阵，因为比较普遍。

很显然，相机的变换可以分成两个阶段：旋转和平移。我们先来看旋转。

在空间坐标系中，相机旋转之前世界坐标系xyz和相机坐标系u0v0n0的各个轴向量的方向相同，有关系： | u0 | | x |P = |Pu0 Pv0 Pn0| x | v0 | = |Px Py Pz| x | y |

| n0 | | z |

相机和向量P都旋转之后，有关系： | u | | x |P' = |Pu0 Pv0 Pn0| x | v | = |Px' Py' Pz'| x | y |

| n | | z |

现在看

| u | | x ||Pu0 Pv0 Pn0| x | v | = |Px' Py' Pz'| x | y |

| n | | z |

因为|x y z|^T为一个单位阵，且Pu0 = Px， Pv0 =Py， Pn0 = Pz。所以得到

| u | |Px Py Pz| x | v | = |Px' Py' Pz'|

| n |

即|Px Py Pz|和相机一起旋转后变成|Px' Py' Pz'|，即P x R = P'，而旋转变换矩阵R就是：

| u |

| v |

| n |

写成标准4x4矩阵：

| ux uy uz 0|

| vx vy vz 0|

| nx ny nz 0|

| 0 0 0 1|

平移矩阵T很简单：

| 1 0 0 0 |

| 0 1 0 0 |

| 0 0 1 0 |

| x y z 1 |

则相机矩阵就是：

| ux uy uz 0 | | 1 0 0 0 | | vx vy vz 0 | | 0 1 0 0 |C = R x T = x | nx ny nz 0 | | 0 0 1 0 | | 0 0 0 1 | | x y z 1 |

它的逆矩阵，即相机的逆变换矩阵为