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矩阵

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

相机

渲染需要定义相机参数，以确定观察位置和方向。最基本的参数是相机位置与视线方向。以下示例定义了相机的位置和观察目标：

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 方向相反，可由 eye 指向目标的向量反向得到：

// 归一化得到单位方向向量，便于后续计算。
const zAxis = normalize(subtract(eye, target));

project

（`eye` 方向向量减去 `lookAt` 方向向量；图中为归一化前的结果）

已知相机坐标系的 +z 轴后，还需要指定世界坐标系的up向量（表示世界空间的向上方向），通常取（0，1，0）。需注意该up向量并非相机坐标系的 +y 轴；后文将根据已知信息推导计算相机坐标系的 +y 轴。

根据右手坐标系规则，相机坐标系的 +x 轴可由世界坐标系的 up向量与相机坐标系的 +z 轴叉乘得到：

特殊情况下当 世界坐标系的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 轴： 
const yAxis = normalize(cross(zAxis, xAxis));

project

（+y 轴归一化前的结果）

得到相机坐标系的三个轴后，需要将顶点坐标转换到该坐标系。向量的各分量即其在相应坐标轴上的投影长度，因此只需计算顶点在相机坐标系各轴上的投影。

点乘公式为：

a⋅b=∥a∥∥b∥cos⁡θa \cdot b = \lVert a \rVert \lVert b \rVert \cos \theta

其中，a 为顶点向量，b 为相机坐标系的轴向量。轴向量已归一化，模长为 1，因此：

a⋅b=∥a∥cos⁡θa \cdot b = \lVert a \rVert \cos \theta

该结果即为顶点在对应轴向量上的投影值。

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 轴方向上的近平面。位于近平面之前的物体不会进入成像范围。

远平面zFar

成像范围在 z 轴方向上的远平面。位于远平面之后的物体不会进入成像范围。

成像宽高比aspectRatio

成像视野的宽高比例。

视野角度eyeFov

成像视野上下边界之间的夹角。

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

将透视锥体挤压成长方体

project 透视投影的成像视野如上图所示。将近平面与远平面之间的视野范围压缩为下图红色的长方体区域。 project

该压缩操作具有以下特性：

project

（透视投影的侧面图）

压缩操作对应一个 4×4 矩阵，输入和输出均为（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。因此：

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’）转换为点向量。为使 w’ 取正值，可在投影前将顶点乘以下列矩阵，对 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（zFar 为正数）时，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; // 近平面的下边界

为将成像区域中心平移至原点，先根据已知条件求出成像范围中心：

(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 表示成像范围的水平宽度；该长方体最终缩放为边长为 2 的立方体。对应的缩放矩阵为：

[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}