OpenGL坐标变换及其数学原理，两种摄像机交互模型（附源程序）

实验平台：win7，VS2010

先上结果截图（文章最后下载程序，解压后直接运行BIN文件夹下的EXE程序）：

a.鼠标拖拽旋转物体，类似于OGRE中的“OgreBites::CameraStyle::CS_ORBIT”。

      

b.键盘WSAD键移动镜头，鼠标拖拽改变镜头方向，类似于OGRE中的“OgreBites::CameraStyle::CS_FREELOOK”。

      

1.坐标变换的一个例子，两种思路理解多个变换的叠加

现在考虑Scale(1,2,1); Transtale(2,1,0); Rotate(pi/4,(0,0,1)); 这3个变换（下文用S, T, R简写），作用到原先中心位于原点边长为2的立方体上的情况。

坐标系显示说明及变换前的场景如下：

以上变换用OpenGL（经典管线）和GLM实现代码分别如下：

glMatrixMode(GL_MODELVIEW);
glPushMatrix();
    glScalef(1, 2, 1);
    glTranslatef(2, 1, 0);
    glRotatef(45, 0, 0, 1);
    glutSolidCube(2);
    draw_frame(1.5f);
glPopMatrix();

glm::mat4 t = glm::scale( glm::vec3(1,2,1) )
     * glm::translate( glm::vec3(2,1,0) )
     * glm::rotate( 45.0f, glm::vec3(0,0,1) );
glMatrixMode(GL_MODELVIEW);
glPushMatrix();
    glMultMatrixf(&t[0][0]);
    glutSolidCube(2);
    draw_frame(1.5f);
glPopMatrix();

变换后的场景如下图：

现在，可以用两种思路来理解S(1,2,1); T(2,1,0); R(pi/4,(0,0,1)); 这三个变换的叠加。

全局坐标变换，所有的变换在一个全局的固定的坐标系下进行，所有操作均以这个坐标系为参考，注意缩放相对于原点进行，这时的变换顺序和代码顺序正好相反，为 R(pi/4,(0,0,1)); T(2,1,0); S(1,2,1); ；

坐标系变换，变换针对坐标系框架进行，所有操作以当前坐标系为参考，所有变换施加后得到一个新坐标系，在这个坐标系中绘制物体，这时的变换顺序和代码相同，为 S(1,2,1); T(2,1,0); R(pi/4,(0,0,1)); ；

两种的图示如下：

这里有几点需要说明或者强调一下。思路1全局坐标变换（左图），最后一步S(1,2,1);相对于全局坐标系的原点，而不是物体中心，所以物体的中心发生了变化。思路2物体坐标系变换（右图），第2步T(2,1,0);以物体坐标系为参考，因为物体坐标系的Y轴在上次变换中被拉长了，所以Y轴的1长度也被拉长了，第3步R(pi/4,(0,0,1));也以物体坐标系为参考，因为物体坐标系的Y轴被拉长了，顶点的旋转轨迹在我们看来是个椭圆，如图中青色所示，同样，所旋转的45度在我们看来也不是45度（物体坐标系并不知道这一点，它只根据当前坐标系进行变换，并且总是觉得自己的旋转轨迹是圆，角度是45度）。

这两种思路显然不是巧合，它们的背后有深刻的数学原理，请接着看下一节。

2.坐标变换的数学原理，两种思路背后的数学解释

因为涉及好多数学公式，这里采用一种新的撰文形式，即PPT加讲解的形式，每页PPT截图后面有旁白解释。如果嫌截图不清楚，文章最后给出了PPT的下载链接。

数域一般取实数集或复数集。基是一组线性无关的向量组，而不一定是互相正交（垂直）的。坐标用列向量表示，在OpenGL中点的坐标也用列向量表示。

基的定义，空间中任一向量均能用基线性表示，另一组基中的向量也如此。T为n介方阵。注意这里T是乘在右边。

Y的公式同理直接写出来了。X, Y均为列向量。注意这里T乘在Y的左边，因为Y是列向量嘛，对比前一页PPT，基变换公式中T乘在左边，这种差异是关键，且往下看。目前讲到的基变换与坐标变换均为线性变换（符合f(ax+by)=af(x)+bf(y)的称为线性），线性变换将原点变换为原点，而OpenGL中的变换可以平移，下面讲到，这是仿射变换。

R表示实数集合，|A|表示A的行列式（determinant，有时也表示为det(A)）。限制A的行列式不为0是要求A非奇异（not singular, invertible，可逆），因为A不可逆时可能将直线映射为一点，即将n维空间压缩为小于n维。另外A的行列式如果为负则变换产生镜像（如将右手系变换为左手系）。之前用X,Y表示坐标也即列向量，现在用粗体小写字母x,y,b表示列向量。PPT中用到了分块矩阵表示，注意A为n×n，x,y,b为n×1，粗体0是1×n个0。x^T表矩阵转置（transposition）。扩充第n+1个坐标是为了能够表示平移，也就是说n+1维空间的线性变换可以表示n维空间的平移。第n+1个坐标还可以用来分辨n维空间中的点与向量（或者说是方向），即第n+1个坐标不为0时表示点，为0时表示向量，不为0且不为1时要将所有坐标都缩放一个倍数使之为1。第n+1个坐标为0时可以从两个角度理解，一是理解成两个点的差，点的第n+1个坐标都是1，做差后第n+1个坐标为0，两个点的差也就是向量，二是将其理解为第n+1个坐标w是从1逼近0，这时可以表示无穷远处的点，也就是一个方向。注意这里的变换矩阵T有固定的形式，即最后一行为n个0接1个1，仿射变换只是n+1维空间的特殊线性变换（自由度小于(n+1)²小于等于n²+n）。如果T的最后一行的前n个元素不为0，那么变换可能将直线变为曲线（请自行举例），即变换后的坐标是原坐标的有理分式（这在OpenGL投影矩阵中被应用）。

原基为向量，前n个元素是齐次坐标系中的向量，这里将线性变换（没有平移）推广到仿射变换，即加入原点，原点是新基中唯一的点（其他为向量）。

这里将之前用的字母T改用A，现在的T表示齐次坐标下的变换矩阵（见PPT第2页）。注意b其实是第一组基下的坐标(和A一样)，这组坐标和基相乘得到它表示的向量。这里再次注意T在基变换和坐标变换公式中的位置。强调一下，T并不是自由的n+1介方阵，它的最后一行固定为n个0接1个1。第n+1个坐标w为0时表示的向量（认为是两个点的差），向量的仿射变换可以看成是其两个端点仿射变换后做差，这时T的平移部分将被抵消（请看下一节仿射变换的分解），也就是说w为0的向量的仿射变换只和T的旋转和缩放部分有关（也就是自由向量的概念，向量只有方向和大小，没有起点）。

这里顺便提一下，矩阵相乘的几何意义就是变换的叠加，即线性映射的叠加。注意一个细节，这里的每个T既表示仿射变换本身，又表示仿射变换的变换矩阵，并没有加以区分，这是合理的，因为仿射变换和仿射变换矩阵之间有一一对应的关系（所有仿射变换构成的空间和所有仿射变换矩阵构成的空间同构）。至此彻底了解了两种思路的数学原理。再次强调这里的仿射变换T可能不一定是刚体变换，它有可能产生缩放、错切变形。第1节的例子就不是刚体变换，以上的两种思路和解释是对仿射变换成立的，不限于刚体变换（旋转和平移或其叠加）。

3.更深入的数学，坐标变换的分解（矩阵的分解）

接着用PPT的形式~

这里都讲的是三维空间。I表示单位矩阵（identity matrix，数学书中一般用E表示）。||v||表示范数，在向量空间中也就是向量的模长。注意到 (u·x)u=(uu^T)x，u×x（叉乘）等于 u波浪线 矩阵乘 x，旋转公式只要选定 u, u×x, x-(u·x)u 三个新基就很好看懂了（文献[1]第11页）。这里既用字母T表示仿射变换矩阵，又用其表示平移矩阵，T的具体含义可以根据上下文区分不会混淆，用C++术语来说，它们的参数列表不同。旋转矩阵沿xyz轴的特殊形式请自行将v设为特殊值进行推导。那现在的问题是，任意给一个仿射变换矩阵T（要符合最后一行是n个0接1个1），T能否分解为T(x,y,z), S(x,y,z), R(a,(x,y,z))的组合（连乘积）呢？答案是肯定的，请继续往下看。

再次，既用T, R, S表示矩阵，又用其表示平移、旋转、缩放函数，请很据上下文区分。行列式为正的正交矩阵是一个旋转矩阵，对称矩阵是个缩放矩阵（缩放值可能有负值，这时产生镜像，即手性变化），经过对角化后分解为旋转矩阵和沿xyz轴缩放的矩阵（即对角阵）。注意极式分解具有唯一性，对角化不具有唯一性，但不唯一性也仅限于调换对角阵的行或列（相应调换对角阵两边的旋转矩阵）。可以根据T, S, R（平移、缩放、旋转）的逆来构造整个变换T的逆（(AB)^-1=B^-1A^-1，当然也可以不分解直接求逆矩阵）。

4.图形学中的变换模型以及OpenGL的实现

这里讲的变换模型是指一种“思维模型”，也就是说用这个模型去思考可以很方便对物体位置和定向进行操作，而具体的实现能够保证按照这个模型思考一定能够得到正确答案，但这个实现可能根本就不是按部就班的按照模型实现具体坐标的计算，所以还要讲OpenGL的实现。

图形学中的变换模型一般涉及物体坐标系（model space）、世界坐标系（world space）、视觉坐标系（eye space）、规范化设备坐标系（normalized device space）、窗口像素坐标系（window space），这些坐标系中的坐标相应叫做某某坐标，如世界坐标系中的坐标叫做世界坐标（world coordinates）。一个示意图如下（用Blender软件制作和渲染的）：

如图中所标注的，猴头上面的坐标框架表示物体坐标系；那个最大的坐标框架是世界坐标系，水红色的是地板；黑色的摄像机上的是视觉坐标系，视觉坐标系的定义是，镜头所指方向为z负方向，摄像机正上为y正方向，右手法则确定x方向。

坐标的变换如下（请见文献[8]第66页）：

         

1.模型变换，视图变换

现在举例子说明物体坐标到视觉坐标的变换，场景是（请看上面猴头那个图），猴头的中心位于物体坐标系原点，猴头中心位于世界坐标系的(1,1,1)处，摄像机位于世界坐标系的(0,1,5)，摄像机的向上方向沿世界坐标系y正方向，摄像机镜头对准世界坐标系z负方向。对猴头中心来说，它在物体坐标系中坐标(0,0,0)，世界坐标系中坐标(1,1,1)，视觉坐标系中坐标(1,0,-4)。模型变换矩阵为T(1,1,1)，视图变换矩阵为T(0,-1,-5)，如果把模型和视图矩阵合起来就是T(0,-1,-5)T(1,1,1)=T(1,0,-4)（还记得，T(x,y,z)表示平移）。GLM和OpenGL函数中的LookAt函数返回的变换矩阵是，将摄像机设置为函数参数指定的情况所需要的视图矩阵。

2.投影变换

投影变换请见下图（摘自文献[4]）：

          

投影变换后进入坐标裁剪，即落在红色方框外的部分将被裁剪掉。

3.透视除法

投影变换后齐次坐标的第4个分量w可能不为1，透视除法即将xyz分量都除以w，得到规范化设备坐标（特点是xyz分量范围在-1到+1之间），对透视投影而言，这一步是非线性的（远处物体被压缩）。如下图（摘自文献[4]）：

          

   

投影变换和透视除法合起来的效果是，将指定的视景体（也叫平截头体，也就是那个裁剪框）变换为边平行于xyz轴且xyz范围都是-1到+1中心位于原点的正方体。注意z坐标的符号变化。如下图（摘自文献[7]）：

  

4.视口变换

再经过视口变换，即调用OpenGL的glViewport函数，对应到窗口像素，注意，在OpenGL中，像素坐标系的原点位于左下角，向右为x轴正向上为y轴正（而一般图片像素都是以左上角为原点）。具体来说，视口变换将规范化设备坐标的位于[-1,1]之间的z坐标对应到深度值，一般在[0,1]（值越小离摄像机越近，z=-1对应d=0，z=+1对应d=1，d为深度值），将(-1,-1,z)对应到屏幕(0,0,d)点，其中d为深度值，将(1,1,z)对应到屏幕(w,h,d)点，其中w,h为窗口的宽和高，其他点按线性插值。如下图（摘自文献[4]）：

5.OpenGL实现

具体到OpenGL的实现，OpenGL和数学中相同采用右手系，OpenGL把模型变换和视图变换合二为一，即模型视图矩阵。OpenGL和GLM的变换矩阵都是按照列优先存储在内存中，这和C++二维数组不同，其实，GLM中的4x4矩阵是由4个列向量组成的。按照上面的分析，当OpenGL的模型视图矩阵和投影矩阵均为单位阵时，这时摄像机位于世界坐标系原点看向z负方向，向右方向沿x轴正方向，向上方向沿y正方向，由于投影矩阵为单位阵，这时为正交投影（另一种是透视投影），裁剪面为xyz的±1，也就是说，对应到最后的显示窗口，x方向向右，y方向向上，z方向垂直屏幕向外，窗口中心对应坐标原点，窗口边缘对应±1，并且z值小的片断遮挡z值大的片断（正好和离摄像机的远近关系反了，这是因为没有对z坐标进行变号）。对了，OpenGL除了模型视图矩阵和投影矩阵之外，还有文理坐标变换矩阵和颜色变换矩阵。请见OpenGL官方手册文献[8]2.12和2.16。

5.两种摄像机交互模型

现在用前面的知识实现两种最常见的摄像机交互模型，先说下对上面说的变换模型的实现，程序有如下全局变量：

glm::mat4 transform_camera(1.0f); // 摄像机的位置和定向，即摄像机在世界坐标系中位置
glm::mat4 transform_model(1.0f);  // 模型变换矩阵，即物体坐标到世界坐标
glm::vec4 position_light0(0);     // 光源位置，世界坐标系中的坐标
float speed_scale=0.1f;           // 鼠标交互，移动速度缩放值

在绘制函数中，这些全局变量被应用如下（第一行之所以求逆，是因为model_view_matrix表示的是视觉坐标到世界坐标的变换矩阵，也就是摄像机在世界坐标系中的位置，这里需要的是将世界坐标变换到视觉坐标）：

glm::mat4 model_view_matrix = glm::affineInverse(transform_camera);
glMatrixMode(GL_MODELVIEW);
glLoadMatrixf(&model_view_matrix[0][0]);
glLightfv(GL_LIGHT0, GL_POSITION, &position_light0[0]); // 位置式光源
draw_world(10,3, true, true, true); // 绘制世界

model_view_matrix *= transform_model;
glMatrixMode(GL_MODELVIEW);
glLoadMatrixf(&model_view_matrix[0][0]);
draw(); // 绘制物体

第一种拖拽球模型（姑且叫做拖拽球吧），它设想窗口中心有个虚拟的球，假设鼠标位于这个球面靠近我们的半面，即z≥0的半球面（使用MeshLab软件制作）：

数学描述如下，图中的 坐标系是视觉坐标系，对应到屏幕也就是向右为x轴，向上为y轴，垂直屏幕向外为z轴：

推导旋转矩阵如下：

这里都假设 a, b 两点在距离屏幕中心小于r的情况，如果大于r（其实是大于水红色圆半径）请见上图的中的 p 点，用 p撇 代替 p 点，水红色圆半径小于r是希望鼠标在距离中心大于r的地方沿径向移动物体也能产生旋转。OpenGL官方Wiki上有个更好的解决方法，见文献[9]。

第二种漫游模型（姑且叫漫游吧），要求按下键盘WSAD键摄像机前进、后退、左移、右移，鼠标左右移动时摄像机镜头左右扫动，鼠标上下移动时摄像机镜头做俯仰动，如下图所示：

注意鼠标左右移动的旋转轴要沿世界坐标系的y轴旋转，而不是摄像机自己的y轴，以防止视角倾斜，鼠标上下移动就沿摄像机自己的x轴旋转就行。键盘WSAD键沿摄像机的z轴和x轴移动就行。

下面是程序实现，程序中的键盘和鼠标响应如下：

1.WS键，摄像机沿视觉坐标系z轴移动，AD键，摄像机沿视觉坐标系x轴移动；

transform_camera *= glm::translate( speed_scale*glm::vec3(dx,0,dz) );

2.上下键，摄像机沿世界坐标系y轴移动，这里v为世界坐标系的y轴单位向量（不是点，所以第四个分量为0，这点很重要，若写成vec4(0,1,0,1)将得到错误结果）在视觉坐标系中的坐标；

glm::vec3 v = glm::vec3( glm::affineInverse(transform_camera)*glm::vec4(0,1,0,0) );
transform_camera *= glm::translate( speed_scale * dy * v );

3.左右键，摄像机沿视觉坐标系z轴旋转；

transform_camera *= glm::rotate( speed_scale*dx, glm::vec3(0,0,1) );

4.鼠标右键拖拽，上下移动时摄像机沿视觉坐标x轴旋转，左右移动时摄像机沿世界坐标系y轴转动；

transform_camera *= glm::rotate( speed_scale*dy, glm::vec3(1,0,0) );
glm::vec3 v = glm::vec3( glm::affineInverse(transform_camera)*glm::vec4(0,1,0,0) );
transform_camera *= glm::rotate( -speed_scale*dx, v );

5.鼠标左键拖拽，物体按拖拽球旋转；

void drag_ball(int x1, int y1, int x2, int y2, glm::mat4& Tmodel, glm::mat4& Tcamera)
{
    float r = (float)std::min(win_h, win_w)/3;
    float r2 = r*0.9f;
    float ax = x1-(float)win_w/2, ay = y1-(float)win_h/2;
    float bx = x2-(float)win_w/2, by = y2-(float)win_h/2;
    float da = std::sqrt(ax*ax+ay*ay), db = std::sqrt(bx*bx+by*by);
    if(std::max(da,db)>r2){
        float dx, dy;
        if(da>db){ dx = (r2/da-1)*ax; dy = (r2/da-1)*ay;
        }else{     dx = (r2/db-1)*bx; dy = (r2/db-1)*by; }
        ax += dx; ay +=dy; bx += dx; by += dy;
    }
    float az = std::sqrt( r*r-(ax*ax+ay*ay) );
    float bz = std::sqrt( r*r-(bx*bx+by*by) );
    glm::vec3 a = glm::vec3(ax,ay,az), b = glm::vec3(bx,by,bz);
    float theta = std::acos(glm::dot(a,b)/(r*r));
    glm::vec3 v2 = glm::cross(a,b);
    // v2是视觉坐标系中的向量，v是v2在物体坐标系中的坐标
    glm::vec3 v = glm::vec3(
        glm::affineInverse(Tmodel) * Tcamera * glm::vec4(v2[0],v2[1],v2[2],0) );
    Tmodel *= glm::rotate( theta*180/3.14f, v );
}

6.鼠标中键拖拽，相当于AD键和上下键；

7.鼠标中键滚动，相当于WS键。

以上代码，以可读性和方便说明原理为目标，所以实现上不很高效，尤其是用transform_camera表示摄像机位置和定向而不是视图矩阵，导致每次都要求transform_camera的逆，可以利用(AB)^-1=B^-1A^-1等公式进行等价变换提高效率。

6.进阶，变换的插值

很多时候，我们希望对变换进行插值，比如，指定物体在开始和结束两个时刻的位置和定向（即物体的transformation），希望在这两个时间点的中间时刻物体能够平滑的变换，从而实现关键帧动画，再比如，我们指定开始和结束两个时刻的摄像机的transformation，希望摄像机的transformation能够被插值，从而实现视角的平滑变化。这个问题可以归结为T(0)=T_begin, T(1)=T_end，求T(t), 0<t<1，使得T(t)随着t平滑变化，这个问题并不像想象中那么简单，T(t)=tT_begin+(1-t)T_end这个函数并不能做到定向（旋转）的平滑变化，甚至都做不到保持物体形状不变（刚体变换）。解决方法涉及高深的数学知识，如矩阵的指数和对数，甚至是群论和李代数，请参考文献[1]。

源程序下载：链接http://pan.baidu.com/s/1hqrG98K 密码: jmc5

PPT下载（如果下载后显示要修复，请右键文件，属性，点下面解除锁定按钮）：链接http://pan.baidu.com/s/1c0lJigw 密码: isds

参考文献

1. Ochiai, H. and Anjyo, K., Mathematical basics of motion and deformation in computer graphics. in ACM SIGGRAPH 2014 Courses, (Vancouver, Canada, 2014), ACM, 1-47.（到ACM网站下载，可能需要大学IP）；
2. 《高等代数简明教程》（上册，第二版，蓝以中编著，北京大学出版社，2007），第4章（到当当网买）；
3. Angle E, Shreiner D. Interactive Computer Graphics: A Top—down Approach with Shader-based OpenGL. 2011. Chapter 3（到亚马逊买）；
4. http://www.opengl-tutorial.org/, Tutorial 3.
5. OpenGL Mathematics (GLM), Open source library.
6. 《OpenGL编程指南》（原书第7版，Dave Shreiner等著，李军等译，机械工业出版社，2011），第3章（到当当网买）；
7. University of Freiburg的Computer Graphics小组课程主页（去课程主页，图形管线PPT）;
8. OpenGL 3.3 Compatibility Profile Specification (updated March 11, 2010), 2.12, 2.16（去官网下载）;
9. https://www.opengl.org/wiki/Trackball。