实现MVP矩阵

在图形渲染管线中，MVP矩阵的功能将模型自身坐标系一步步变换到最终可以显示在屏幕上的坐标的过程，本篇文章将深入分析MVP矩阵的实现。

阅读本篇文章需要有使用gl-matrix等类似第三方矩阵库进行3d渲染的前置知识 https://developer.mozilla.org/en-US/docs/Web/API/WebGL_API/Tutorial

本篇文章配套例程序：https://github.com/xxkl1/mvp_matrix

MVP矩阵

MVP矩阵分别由Model Matrix、View Matrix、和Projection Matrix组成

Model Matrix: 将模型放置在世界空间中的期望的位置；

View Matrix: 世界空间转为相机空间；

Projection Matrix： 实现正交投影或者透视投影，并转换为标准屏幕空间；

Model Matrix

世界空间

Model Matrix前需要了解世界空间，世界空间和我们身处的世界相似，其范围是无限大，即x，y，z的正负轴取值可以无限大和无限小，并且坐标系的原点在（0, 0, 0);

Model Matrix的作用

Model Matrix是将模型通过变换，放置到我们想要的位置，如果想要放大缩小等，还需要涉及缩放等线性操作操作。

下面是一个正方体模型的顶点数据，其定义了正方体的6个面，每个面有4个顶点，顶点使用（x，y，z）表示

const positions =  [
    // 前面
    -3, -3, 3,
    3, -3, 3,
    3, 3, 3,
    -3, 3, 3,

    // 后面
    -3, -3, -3,
    -3, 3, -3,
    3, 3, -3,
    3, -3, -3,

    // 顶面
    -3, 3, -3,
    -3, 3, 3,
    3, 3, 3,
    3, 3, -3,

    // 底面
    -3, -3, -3,
    3, -3, -3,
    3, -3, 3,
    -3, -3, 3,

    // 右面
    3, -3, -3,
    3, 3, -3,
    3, 3, 3,
    3, -3, 3,

    // 左面
    -3, -3, -3,
    -3, -3, 3,
    -3, 3, 3,
    -3, 3, -3,
]

model

（模型的顶点数据，模型原点一般是（0，0，0））
有了模型数据后，需要将该模型放置在世界空间期望的位置，类似我们需要将一个物体放置在房间的某个地方一样

translate

（调整模型放置在世界空间的位置）

不仅是位置，还期望给模型做一些旋转操作，类似希望模型的哪个面贴在房间墙面 model

（调整实现模型围绕x轴，y轴，z轴旋转）

除了平移，旋转，还有缩放等操作，这些其实都是在世界空间下，对模型做线性变换操作，使得模型转换成最终在世界空间中的期望的状态。

View Matrix

相机

在世界中，我们需要相机/人眼来看这个世界，以实现成像，所以需要定义相机的相关参数，相机的参数主要有两个，相机位置、相机的视线方向，例如下面就定义了相机在（0，6，6）的位置，看向（0，0，1）

const eye = {
  position: [1, -1, -1],
  lookAt: [2, -2, -2],
}

相机空间

在进行相机成像前，需要将世界空间转移到相机空间，以方便后面的计算。相机空间的特征是：

  1. 相机的视线方向归一化后是（0， 0， -1），即视线方向是-z轴；
  2. 相机位置在坐标系原点（0， 0， 0）；

世界空间转相机空间思路

  1. 改变相对坐标系 先找到在世界空间中，相机坐标系+x、+y和+z轴的方向 由于相机坐标系中的视线方向的-z轴方向，所以+z轴的方向和lookAt的方向相反，即得到相机的+z轴是eye方向向量减去lookAt方向的向量
// normalize是为了得到单位向量，因为结果仅代表方向，也是方便后面的运算
const zAxis = normalize(subtract(eye, target));

project

（eye方向向量减去lookAt方向的向量，归一化前的结果如上图所示）

现在知道了相机坐标系的-z轴，还有世界坐标系的+y轴，即up向量（这个一般是（0， 1， 0）），根据叉乘定义，up和z轴的叉乘，可以得到相机坐标系的+x轴

（根据叉乘+右手坐标系定律，右手第一个手指拇指指向up，第二个手指食指指向相机坐标+z轴，那么第三个手指中指的方向就是相机坐标系的+x轴，且相机+x轴和相机坐标系的+z轴、up都垂直）

这里还有一个细节，当up和相机+z轴平行的时候，会出现叉乘结果运算错误的情况，这个时候需要进行up方向修正（使用（0，0，1））并进行叉乘的重新计算

let up = new Vec4([0, 1, 0, 1]);
let xAxis = normalize(cross(up, zAxis));
if (isNaN(xAxis.value[0]) || isNaN(xAxis.value[1]) || isNaN(xAxis.value[2])) {
    up = new Vec4([0, 0, 1, 1]);
    xAxis = normalize(cross(up, zAxis));
}

project

（+x轴归一化前的结果如上图所示）
知道了相机坐标系的+z轴和+x轴的方向，同样根据叉乘+右手坐标系定律，可以得到相机坐标系的+y轴

（右手第一个手指拇指指向+z轴，第二个手指食指指向相机坐标+x轴，那么第三个手指中指的方向就是相机坐标系的+y轴）

const yAxis = normalize(cross(zAxis, xAxis));

project

（+y轴归一化前的结果如上图所示）

得到了相机坐标系的三个轴后，需要将顶点坐标的相对坐标系转换到相机坐标系上，已知向量的三维分量其实就是点在三维轴的投影点到原点的距离，因此需要求的值就变成，计算顶点在相机坐标系下的投影

点乘公式如下

a⋅b=∥a∥∥b∥cosθa⋅b=∥a∥∥b∥cosθ

a是顶点的向量值，b是相机坐标系的轴向量，轴向量的模进行了归一化，所以都是1，得到

a⋅b=∥a∥cosθa⋅b=∥a∥cosθ

刚好就是顶点在该轴向量的投影值

project

（顶点在轴向量上的投影）

在笛卡尔坐标系下，向量的点乘公式是

a⃗⋅b⃗=(xayaza)⋅(xbybzb)=xaxb+yayb+zazb\vec{a}\cdot\vec{b}= \begin{pmatrix}x_a\\y_a\\z_a\end{pmatrix}\cdot \begin{pmatrix}x_b\\y_b\\z_b\end{pmatrix} = x_a x_b + y_a y_b + z_a z_b

得到顶点在x轴向量的投影是

xAxis⃗⋅point⃗=(xAxis.xxAxis.yxAxis.z)⋅(point.xpoint.ypoint.z)=xAxis.x∗point.x+xAxis.y∗point.y+xAxis.z∗point.z\vec{xAxis}\cdot\vec{point}= \begin{pmatrix}xAxis.x\\xAxis.y\\xAxis.z\end{pmatrix}\cdot \begin{pmatrix}point.x\\point.y\\point.z\end{pmatrix} = xAxis.x * point.x + xAxis.y * point.y + xAxis.z * point.z

y轴向量和z轴向量以此类推

最终的x轴向量、y轴向量，z轴向量投影矩阵表示就是

∣xAxis.xxAxis.yxAxis.z0yAxis.xyAxis.yyAxis.z0zAxis.xzAxis.yzAxis.z00001∣∣point.xpoint.ypoint.z1∣\begin{vmatrix} xAxis.x & xAxis.y & xAxis.z & 0 \\ yAxis.x & yAxis.y & yAxis.z & 0 \\ zAxis.x & zAxis.y & zAxis.z & 0 \\ 0 & 0 & 0 & 1 \end{vmatrix} \begin{vmatrix} point.x \\ point.y \\ point.z \\ 1 \end{vmatrix}

在相机空间，相机位置在原点，所以顶点需要进行整体的平移变换

∣100−eye.x010−eye.y001−eye.z0001∣\begin{vmatrix} 1 & 0 & 0 & -eye.x \\ 0 & 1 & 0 & -eye.y \\ 0 & 0 & 1 & -eye.z \\ 0 & 0 & 0 & 1 \end{vmatrix}

注意不能进行相对坐标系变换，再平移，会因为前后坐标系的不同，出现偏移的情况 （但是为什么在计算相机坐标系轴向量时候，不需要先进行平移变换？原因是轴向量是一个方向，在同样的世界空间下，平移前平移后，方向是不变的，所以计算相机坐标系轴向量的时候不需要平移））

最终的View Matrix就是

∣xAxis.xxAxis.yxAxis.z0yAxis.xyAxis.yyAxis.z0zAxis.xzAxis.yzAxis.z00001∣∣100−eye.x010−eye.y001−eye.z0001∣=∣xAxis.xxAxis.yxAxis.z−dot(eye,xAxis)yAxis.xyAxis.yyAxis.z−dot(eye,yAxis)zAxis.xzAxis.yzAxis.z−dot(eye,yAxis)0001∣\begin{vmatrix} xAxis.x & xAxis.y & xAxis.z & 0 \\ yAxis.x & yAxis.y & yAxis.z & 0 \\ zAxis.x & zAxis.y & zAxis.z & 0 \\ 0 & 0 & 0 & 1 \end{vmatrix} \begin{vmatrix} 1 & 0 & 0 & -eye.x \\ 0 & 1 & 0 & -eye.y \\ 0 & 0 & 1 & -eye.z \\ 0 & 0 & 0 & 1 \end{vmatrix} = \begin{vmatrix} xAxis.x & xAxis.y & xAxis.z & -dot(eye, xAxis) \\ yAxis.x & yAxis.y & yAxis.z & -dot(eye, yAxis) \\ zAxis.x & zAxis.y & zAxis.z & -dot(eye, yAxis) \\ 0 & 0 & 0 & 1 \end{vmatrix}

Project Matrix

透视投影和正交投影

project

（左：透视投影，右：正交投影，图片来源：GAMES101）

和透视投影相对的是正交投影，正交投影和投影投影最大的区别是透视具有 近大远小 的效果，而正交投影无论多远的物体，成像和原本的物体宽高不会发生变化，所以透视投影会更符合我们真实生活的成像效果。

如何实现近大远小？

那么具体透视投影是如何实现近大远小？在上面投影图中，可以看到透视投影相机到物体的视野范围是一个锥体，视野范围和相机离物体的距离成正比。而正交投影的视野范围是固定的，呈现一个长方体。如果我们将投影中的锥体挤压成一个和正交投影一样的长方体，那么视野范围越大的，被挤压的量就越大，从而物体就越小。 project

（近大远小，图片来源：pixnio）

上图中，靠近相机的视野范围小，远离相机视野范围大，视野范围和被挤压量成正比，远处的山被挤压的量比人大，使得远处山峰最终成像比人还小。

相机位置

开始推导前，先确定相机位置和相机视线中心朝向

相关参数

近平面zNear

进行成像的z轴近平面，小于zNear的物体将无法被纳入成像范围。

远平面zFar

进行成像的z轴远平面，大于zNear的物体将无法被纳入成像范围。

成像宽高比aspectRatio

成像视野的宽高比。

视野角度eyeFov

成像视野范围最高点和最低点的夹角。

透视投影矩阵具体推导过程

将透视锥体挤压成长方体

project 透视投影的成像视野如上图所示，近平面和远平面的视野范围挤压成下面红色的长方体区域 project

该挤压操作有几个特性

project

（透视投影的侧面图）

挤压操作对应的是一个4x4矩阵，输入和输出都是（x, y, z, w) 齐次坐标，a-p都是未知

∣abcdefghijklmnop∣\begin{vmatrix} a & b & c & d \\ e & f & g & h \\ i & j & k & l \\ m & n & o & p \end{vmatrix}

使用特殊点法来进行矩阵推导，特殊点取上图的成像范围内某个面的最上面的点(x, y, z, w)，其挤压后点是(x’, y’, z’, w’)

相机和近平面形成的三角形，和相机到被成像物体面构成的三角形呈相似三角形，根据相似三角形法则，对应的边的比例一致和挤压操作里的已知特性1可以得到

yy′=zzNear\frac{y}{y'} = \frac{z}{zNear}

可以推导出

y′=y⋅zNearzy' =\frac{y \cdot zNear}{z}

同理可以得到

x′=x⋅zNearzx' =\frac{x \cdot zNear}{z}

在齐次坐标下，x，y，z的值最终需要除于w，所以需要将上面的x，y，z除于w分量，得出

y′w′=yw⋅zNearzw\frac{y'}{w'} =\frac{\frac{y}{w} \cdot zNear}{\frac{z}{w}}

和

x′w′=xw⋅zNearzw\frac{x'}{w'} =\frac{\frac{x}{w} \cdot zNear}{\frac{z}{w}}

不妨设挤压前的w分量都是1，那么可以得到

y′w′=y⋅zNearz\frac{y'}{w'} =\frac{y \cdot zNear}{z}

和

x′w′=x⋅zNearz\frac{x'}{w'} =\frac{x \cdot zNear}{z}

设变换后的w’分量是z值，可以得到

y′=y⋅zNeary' = y \cdot zNear

和

x′=x⋅zNearx' = x \cdot zNear

但是实际这里会有点问题，由于相机是在原点，向负z轴看，所以成像范围内物体的z都是负值。而最后(x’, y’, z’, w’)需要通过(x’/w’, y’/w’, z’/w’)转成点向量的。所以，需要在进行投影前，先将顶点的乘上下面矩阵，对z进行取反，使得z变成正数，相机视线也变成+z轴，这样相机的成像也是和原本是一样的，但是得注意，在进行深度测试的时候，需要变成z值越小，越靠近相机的情况。

∣1000010000−110001∣\begin{vmatrix} 1 & 0 & 0 & 0 \\ 0 & 1 & 0 & 0 \\ 0 & 0 & -1 & 1 \\ 0 & 0 & 0 & 1 \\ \end{vmatrix}

根据矩阵乘法求x’可以得到

a⋅x+b⋅y+c⋅z+d⋅w=x′a \cdot x + b \cdot y + c \cdot z + d \cdot w = x'

根据已知道

x′=x⋅zNearx' = x \cdot zNear

可以推导出，a = zNear，b、c和d都等于0

根据矩阵乘法求y’可以得到

e⋅x+f⋅y+g⋅z+h⋅w=y′e \cdot x + f \cdot y + g \cdot z + h \cdot w = y'

根据已知道

y′=y⋅zNeary' = y \cdot zNear

可以推导出，f = zNear，e、g和h都等于0

根据矩阵乘法求z’可以得到

i⋅x+j⋅y+k⋅z+l⋅w=z′i \cdot x + j \cdot y + k \cdot z + l \cdot w = z'

设i和j都等于0，并且w = 1，所以变成

k⋅z+l=z′k \cdot z + l = z'

进行齐次坐标除于w’分量转为点向量

k⋅z+lw′=z′w′\frac{k \cdot z + l}{w'} = \frac{z'}{w'}

由于w’ = z，可以得到

k⋅z+lz=z′w′\frac{k \cdot z + l}{z} = \frac{z'}{w'}

根据 挤压操作的特性3，当z = zNear（zNear是正数）时候，z’/w’等于zNear，即

k⋅zNear+lzNear=zNear\frac{k \cdot zNear + l}{zNear} = zNear

同理，当z = zFar（zNear是正数）时候，z’/w’等于zFar，即

k⋅zFar+lzFar=zFar\frac{k \cdot zFar + l}{zFar} = zFar

化简第一个方程

两边同乘zNear

k⋅zNear+l=zNear2k \cdot zNear + l = zNear^2

得到

l=zNear2−k⋅zNearl = zNear^2 - k \cdot zNear

同理可以得到第二个方向可以简化成

l=zFar2−k⋅zFarl = zFar^2 - k \cdot zFar

可以得到

zNear2−k⋅zNear=zFar2−k⋅zFarzNear^2 - k \cdot zNear = zFar^2 - k \cdot zFar

移项整理得到

k(zFar−zNear)=zFar2−zNear2k(zFar−zNear)=zFar^2−zNear^2

右边进行平方差分解得到

zFar2−zNear2=(zFar−zNear)(zFar+zNear)zFar^2−zNear^2 = (zFar−zNear)(zFar+zNear)

因此得到

k=zFar+zNeark =zFar + zNear

将k值代回l=zNear2−k⋅zNearl = zNear^2 - k \cdot zNear得到

l=zNear2−（zFar+zNear)⋅zNearl = zNear^2 - （zFar + zNear) \cdot zNear

得到

l=−zFar⋅zNearl = -zFar \cdot zNear

根据矩阵乘法求w’可以得到

m⋅x+n⋅y+o⋅z+p⋅w=w′m \cdot x + n \cdot y + o \cdot z + p \cdot w = w'

已知

w′=zw'=z

可以得到o=1，m、n和p等于0

+z轴挤压操作矩阵等于

∣Znear0000Znear0000zFar+zNear−zFar⋅zNear0010∣\begin{vmatrix} Znear & 0 & 0 & 0 \\ 0 & Znear & 0 & 0 \\ 0 & 0 & zFar + zNear & -zFar \cdot zNear \\ 0 & 0 & 1 & 0 \end{vmatrix}

由于前面需要提前进行z轴取负操作，所以最终的挤压矩阵是

∣Znear0000Znear0000zFar+zNearzFar⋅zNear0010∣∗∣1000010000−110001∣=∣Znear0000Znear0000−(zFar+zNear)zFar⋅zNear00−10∣\begin{vmatrix} Znear & 0 & 0 & 0 \\ 0 & Znear & 0 & 0 \\ 0 & 0 & zFar + zNear & zFar \cdot zNear \\ 0 & 0 & 1 & 0 \end{vmatrix} * \begin{vmatrix} 1 & 0 & 0 & 0 \\ 0 & 1 & 0 & 0 \\ 0 & 0 & -1 & 1 \\ 0 & 0 & 0 & 1 \\ \end{vmatrix} = \begin{vmatrix} Znear & 0 & 0 & 0 \\ 0 & Znear & 0 & 0 \\ 0 & 0 & -(zFar + zNear) & zFar \cdot zNear \\ 0 & 0 & -1 & 0 \end{vmatrix}

移动和缩放到标准坐标空间

挤压完成后，需要物体成像的中心点移动坐标原点，和归一化，从而转换为标准坐标空间，以方便进行渲染

const angle = eyeFov / 180.0 * MY_PI; // 角度转换为弧度
const top = zNear * Math.tan(angle / 2); // 近平面上边界
const right = top * aspect_ratio; // 近平面的右边界
const left = -right; // 近平面左边界
const bottom = -top; // 近平面的下边界

要实现物体成像的中心点移动到原点，x轴移动，先根据已知信息求得成像范围中心是

(0,0,zNear+zFar2),(0, 0, \frac{zNear + zFar}{2}),

所以对应平移矩阵是

∣10000100001−(zNear+zFar)20001∣\begin{vmatrix} 1 & 0 & 0 & 0 \\ 0 & 1 & 0 & 0 \\ 0 & 0 & 1 & \frac{-(zNear + zFar)}{2} \\ 0 & 0 & 0 & 1 \end{vmatrix}

将上下左右边界进行缩放归一化得到x，y，z缩放到[-1, 1]，对应缩放矩阵是 r - l是求成像范围水平距离，等等，然后这个长方体需要缩放成222的正方体

[2r−l00002t−b00002zfar−znear00001]\begin{bmatrix} \frac{2}{r - l} & 0 & 0 & 0 \\ 0 & \frac{2}{t - b} & 0 & 0 \\ 0 & 0 & \frac{2}{z_{\text{far}} - z_{\text{near}}} & 0 \\ 0 & 0 & 0 & 1 \end{bmatrix}