Real-Time Rendering 4th Edition学习笔记(三)

Chapter 4: 变换(Transform)

仿射变换(affine transform)是指将一个向量进行一次线性变换后再接一次平移变换。我们用$\mathbf{v}=(v_x\;v_y\;v_z\;0)^T$来表示一个方向向量，用$\mathbf{v}=(v_x\;v_y\;v_z\;1)^T$来表示一个点。

所有的平移、旋转、缩放、反射和剪切矩阵都是仿射的。

基本变换

平移(Translation)

一个点从一个坐标移动向量$\mathbf{t}=(t_x,t_y,t_z)$到另一个坐标可以用如下的平移矩阵来表示：

$$\mathbf{T}(t)=\mathbf{T}(t_{x},t_{y},t_{z})=\begin{bmatrix} 1 & 0 & 0 & t_{x}\\ 0 & 1 & 0 & t_{y}\\ 0 & 0 & 1 & t_{z}\\ 0 & 0 & 0 & 1 \end{bmatrix}$$

图中方块按照平移矩阵$\mathbf{T}(5,2,0)$进行变换，即向x轴正方向移动5个单位，向y正方向移动2个单位

明显可知，一个点$\mathbf{p}=(p_x\;\;p_y\;\;p_z\;\;1)^T$乘以$T(t)$后得到了一个新的点$\mathbf{p{}'}=(p_x+t_x\;\;p_y+t_y\;\;p_z+t_z\;\;1)^T$。一个方向向量$\mathbf{v}=(v_x\;\;v_y\;\;v_z\;\;0)^T$乘以$\mathbf{T}$后保持不变，因为方向向量本身无法进行平移。点和方向都会受平移变换以外其他仿射变换的影响。

$$\mathbf{T}(t)*\mathbf{p}=\begin{bmatrix} 1 & 0 & 0 & t_{x}\\ 0 & 1 & 0 & t_{y}\\ 0 & 0 & 1 & t_{z}\\ 0 & 0 & 0 & 1 \end{bmatrix}*\begin{bmatrix} p_{x}\\ p_{y}\\ p_{z}\\ 1 \end{bmatrix}=\begin{bmatrix} p_{x}+t_{x}\\ p_{y}+t_{y}\\ p_{z}+t_{z}\\ 1 \end{bmatrix}$$

$$\mathbf{T}(t)*\mathbf{v}=\begin{bmatrix} 1 & 0 & 0 & t_{x}\\ 0 & 1 & 0 & t_{y}\\ 0 & 0 & 1 & t_{z}\\ 0 & 0 & 0 & 1 \end{bmatrix}*\begin{bmatrix} v_{x}\\ v_{y}\\ v_{z}\\ 0 \end{bmatrix}=\begin{bmatrix} v_{x}\\ v_{y}\\ v_{z}\\ 0 \end{bmatrix}$$

平移矩阵的逆矩阵为$\mathbf{T}^{-1}(t)=\mathbf{T}(-t)$，即向量t取反。

逆矩阵：给定一个$n*n$的矩阵$\mathbf{A}$，如果存在矩阵$\mathbf{B}$，使$\mathbf{AB}=\mathbf{BA}=\mathbf{I}_n$，则$\mathbf{A}$是可逆的，$\mathbf{B}$为$\mathbf{A}$的逆矩阵，写作$\mathbf{A}^{-1}$。其中$\mathbf{I}_n$为n阶单位矩阵，即形如$\begin{bmatrix}1 & 0 & 0 & ...\\0 & 1 & 0 & ...\\0 & 0 & 1 & ...\\... & ... & ... & 1\end{bmatrix}$的矩阵

旋转(Rotation)

在二维中，旋转矩阵很容易推导。假设现在有一个向量$\mathbf{v}=(v_x,v_y)$，我们将其参数化为$\mathbf{v}=(v_x,v_y)=(rcos\theta,rsin\theta)$。如果我们将这个向量沿逆时针旋转$\phi^{\circ}$，则它将变成$\mathbf{u}=(rcos(\theta+\phi),rsin(\theta+\phi))$。进一步的，我们有：

$$\mathbf{u}=\begin{bmatrix} r\cos(\theta + \phi )\\ r\sin(\theta + \phi ) \end{bmatrix}=\begin{bmatrix} r(\cos\theta \cos\phi -\sin\theta \sin\phi)\\ r(\sin\theta \cos\phi -\cos\theta \sin\phi) \end{bmatrix}\\ =\underbrace{\begin{bmatrix} \cos\phi& -\sin\phi\\ \sin\phi & \cos\phi \end{bmatrix}}_{\mathbf{R}(\phi)}\underbrace{\begin{bmatrix} r\cos\theta\\ r\sin\theta \end{bmatrix}}_{\mathbf{v}}=\mathbf{R}(\phi)\mathbf{v}$$

在三维中，我们通常用$\mathbf{R}_x(\phi)$、$\mathbf{R}_y(\phi)$和$\mathbf{R}_z(\phi)$来分别表示将物体沿x、y、z轴旋转$\phi^{\circ}$。它们的矩阵分别为：

$$\mathbf{R}_{x}(\phi)=\begin{bmatrix} 1 & 0 & 0 & 0\\ 0 & \cos\phi & -\sin\phi & 0\\ 0 & \sin\phi & \cos\phi & 0\\ 0 & 0 & 0 & 1 \end{bmatrix}$$

$$\mathbf{R}_{y}(\phi)=\begin{bmatrix} \cos\phi & 0 & \sin\phi & 0\\ 0 & 1 & 0 & 0\\ -\sin\phi & 0 & \cos\phi & 0\\ 0 & 0 & 0 & 1 \end{bmatrix}$$

$$\mathbf{R}_{z}(\phi)=\begin{bmatrix} \cos\phi & -\sin\phi & 0 & 0\\ \sin\phi & \cos\phi & 0 & 0\\ 0 & 0 & 1 & 0\\ 0 & 0 & 0 & 1 \end{bmatrix}$$

如果将上述4x4的矩阵的最下行和最右列删除，就会得到一个3x3矩阵。这样的3x3矩阵R表示沿着任意轴旋转$\phi^{\circ}$，它的迹(Trace)，即$n*n$矩阵从左上角到右下角对角线上的n个元素之和)为固定值：

$tr(\mathbf{R})=1+2\cos\phi$

$\mathbf{R}_i(\phi)$，即物体围绕i轴旋转$\phi^{\circ}$，其本质就是在i轴上的坐标都保持不变。

所有旋转矩阵的行列式(determinant)都为1，并且都是正交矩阵(orthogonal matrix)。

行列式：二维矩阵$det(\mathbf{A})=\|\mathbf{A}| =det(\begin{bmatrix}a & b\\c & d\end{bmatrix})=ad-bc$，用来表示两个向量所围成的平行四边形的面积。三维则表示三个向量围成的立方体的体积，其余概念类似。

转置矩阵：矩阵$\mathbf{A}$的转置矩阵写作$\mathbf{A}^T$，即将$\mathbf{A}$的行写作$\mathbf{A}^T$的列，$\mathbf{A}$的列写作$\mathbf{A}^T$的行。如： $\begin{bmatrix}1 & 2\\ 3 & 4\\ 5 & 6\end{bmatrix}^T=\begin{bmatrix}1 & 3 & 5\\ 2 & 4 & 6\end{bmatrix}$

正交矩阵：矩阵$\mathbf{Q}$的转置矩阵$\mathbf{Q}^T$和逆矩阵$\mathbf{Q}^{-1}$相同。正交矩阵的行列式必定为1或-1，因为：

$$\because 1=det(I)=det(Q^{-1}Q)=det(Q^TQ)=det(Q^T)*det(Q)=(det(Q))^2$$

$$\therefore det(Q)=\pm 1$$

现在假设我们要将一个物体沿z轴旋转$\phi^{\circ}$，旋转的原点是p。首先，我们要先将物体按照$\mathbf{T}(-p)$进行平移，以保证p点在原点，不受旋转影响。接着，对物体按照$\mathbf{R}_z(\phi)$进行旋转。最后，用$\mathbf{T}(p)$让物体回到原先位置。所以最终的变换矩阵为：

$$\mathbf{X}=\mathbf{T}(\mathbf{p})\mathbf{R}_z(\phi)\mathbf{T}(\mathbf{-p})$$

缩放(Scaling)

缩放矩阵用$\mathbf{S}(s)=\mathbf{S}(s_x,s_y,s_z)$表示，即物体在x、y、z轴方向分别缩放$s_x$、$s_y$、$s_z$倍。

$$\mathbf{S}(s)=\begin{bmatrix} s_x & 0 & 0 &0\\ 0 & s_y & 0 &0\\ 0 & 0 & s_z & 0\\ 0 & 0 & 0 & 1 \end{bmatrix}$$

如果$s_x=s_y=s_z$，则称为均匀(uniform)或各向同性(isotropic)缩放，否则为不均匀(nonuniform)或各向异性(anisotropic)缩放。其逆矩阵为$\mathbf{S}^{-1}(s)=\mathbf{S}(1/s_x,1/s_y,1/s_z)$

当使用齐次坐标时，可以通过修改矩阵最右下角的值来表示缩放。假设我们要将一个物体放大5倍，可以用两种方式来表示：

$$\mathbf{S}=\begin{bmatrix} 5 & 0 & 0 &0\\ 0 & 5 & 0 &0\\ 0 & 0 & 5 & 0\\ 0 & 0 & 0 & 1 \end{bmatrix},\;\;\;\;\;\;\mathbf{S}{}'=\begin{bmatrix} 1 & 0 & 0 &0\\ 0 & 1 & 0 &0\\ 0 & 0 & 1 & 0\\ 0 & 0 & 0 & 1/5 \end{bmatrix}$$

前者可以用来进行不均匀缩放，而后者只能进行均匀缩放。由于后者涉及到齐次坐标的除法，效率可能会降低，除非系统不管右下角的数值是否是1都会进行一次除法计算。

如果s向量中有某个值是负值，则该矩阵为反射矩阵(reflection matrix)或镜像矩阵(mirror matrix)。如果三个参数中有两个是-1，则我们会将物体旋转180°。

反射矩阵通常需要特殊处理。例如，一个三角形三个顶点的顺序是逆时针的，在经过反射矩阵变换后，就变成了顺时针的，如果不处理，就会导致光照和背向剔除出现问题。我们可以通过计算4x4矩阵的左上角3x3区域的矩阵的行列式来判断原矩阵是否是反射矩阵，如果行列式为负数，则说明它是反射矩阵。

剪切(Shearing)

剪切矩阵可以用来扭曲整个场景或模型外观。它一共有6种基本矩阵，分别表示为$\mathbf{H}_{xy}(s),\mathbf{H}_{xz}(s),\mathbf{H}_{yx}(s),\mathbf{H}_{yz}(s),\mathbf{H}_{zx}(s),\mathbf{H}_{zy}(s)$。第一个下标表示在哪个坐标轴上变换，第二个下标表示以哪个坐标轴的值作为系数。剪切矩阵可以表示为：

$$\mathbf{H}_{xz}(s)=\begin{bmatrix} 1 & 0 & s & 0\\ 0 & 1 & 0 & 0\\ 0 & 0 & 1 & 0\\ 0 & 0 & 0 & 1 \end{bmatrix}$$

一个点$\mathbf{p}=(p_x\;\;p_y\;\;p_z\;\;1)^T$乘以上述矩阵后得到了一个新的点$\mathbf{p{}'}=(p_x+sp_z\;\;p_y\;\;p_z\;\;1)^T$

$$\mathbf{H}_{xz}(s)*\mathbf{p}=\begin{bmatrix} 1 & 0 & s &0\\ 0 & 1 & 0 &0\\ 0 & 0 & 1 & 0\\ 0 & 0 & 0 & 1 \end{bmatrix}*\begin{bmatrix} p_x\\ p_y\\ p_z\\ 1\end{bmatrix}=\begin{bmatrix} p_x+sp_z\\ p_y\\ p_z\\ 1\end{bmatrix}$$

其结果如下图所示：

剪切矩阵的逆矩阵为$\mathbf{H}_{ij}^{-1}(s)=\mathbf{H}_{ij}(-s)$

剪切矩阵还有另外一种表现方式是：

$$\mathbf{H}_{xy}{}'(s,t)=\begin{bmatrix} 1 & 0 & s &0\\ 0 & 1 & t &0\\ 0 & 0 & 1 & 0\\ 0 & 0 & 0 & 1 \end{bmatrix}$$

在这种情况下，两个下标的轴向都会发生变换，另外一个轴的值将作为变换的系数

$$\mathbf{H}_{xy}{}'(s,t)*\mathbf{p}=\begin{bmatrix} 1 & 0 & s &0\\ 0 & 1 & t &0\\ 0 & 0 & 1 & 0\\ 0 & 0 & 0 & 1 \end{bmatrix}*\begin{bmatrix} p_x\\ p_y\\ p_z\\ 1\end{bmatrix}=\begin{bmatrix} p_x+sp_z\\ p_y+tp_z\\ p_z\\ 1\end{bmatrix}$$

需要注意的是，所有剪切矩阵的行列式都为1，不会改变物体体积。

变换级联(Concatenation of Transforms)

由于矩阵在乘法运算上的非互换性，变换的级联是和顺序相关的。如下图所示，先旋转再缩放和先缩放再旋转得到的结果完全不同。

将一系列变换矩阵级联成一个矩阵可以提高效率。这个单一矩阵为$\mathbf{C}=\mathbf{TRS}$。这个顺序很重要，一个顶点会先进行缩放，再进行旋转，最后进行位移，即$\mathbf{Cp}=\mathbf{TRSp}=\mathbf{T(R(Sp))}$。

另外需要注意的是，在遵从上述乘法运算顺序的前提下，TRSp可以被自行组合。比如，可以将其组合成(TR)(Sp)。

刚体变换(Rigid-Body Transform)

对于一个坚硬的物体来说，它的变换过程只有位置和角度会发生变化，形状不会发生变化。这种只包含位移和旋转矩阵的变换叫做Rigid-Body Transform，可以表示为：

$$\mathbf{X}=\mathbf{T}(t)\mathbf{R}=\begin{bmatrix} r_{00} & r_{01} & r_{02} & t_x\\ r_{10} & r_{11} & r_{12} & t_y\\ r_{20} & r_{21} & r_{22} & t_z\\ 0 & 0 & 0 & 1\\ \end{bmatrix}$$

它的逆矩阵为$\mathbf{X}^{-1}=(\mathbf{T}(t)\mathbf{R})^{-1}=\mathbf{R}^{-1}\mathbf{T}^{-1}(t)=\mathbf{R}^T\mathbf{T}(-t)$，也即：先将R的左上角3x3矩阵转置，再将T的平移方向取负值，最后用相反顺序将它们俩相乘。

假设我们要和OpenGL Utility Library(GLU)中的_gluLookAt_方法一样来调整相机的朝向。如果相机坐标是$\mathbf{c}$，朝向目标位置在$\mathbf{l}$，相机朝上的方向向量为$\mathbf{u}{}'$，需要求取的是三个基本向量{$\mathbf{r,u,v}$}，如下图所示：

首先计算视角向量(view vector)，$\mathbf{v}=(\mathbf{c}-\mathbf{l})/\| \mathbf{c}-\mathbf{l}\|$，即从目标物体到相机的归一化(normalized)方向向量。向右的方向则为$\mathbf{u}{}'$和$\mathbf{v}$的叉乘的normalize，即$\mathbf{r}=-(\mathbf{v}\times \mathbf{u}{}')/\| \mathbf{v}\times \mathbf{u}{}'\|$。$\mathbf{u}{}'$通常不是正好朝向正上方，所以再做一次叉乘，$\mathbf{u}=\mathbf{v}\times \mathbf{r}$。

在构建相机的变换矩阵$\mathbf{M}$时，我们需要先让相机位移到原点，再使其x、y、z轴分别和$\mathbf{r,u,v}$对齐。如果相机从原点无旋转的状态变换到当前状态经历了T(t)R，那使它的世界坐标转换成视图坐标的矩阵就是：

$$\mathbf{M}=(\mathbf{T}(t)\mathbf{R})^{-1}=\mathbf{R}^{-1}\mathbf{T}^{-1}(t)$$

我们知道，一个矩阵$\begin{bmatrix}x_x & y_x & z_x\\x_y & y_y & z_y\\x_z & y_z & z_z\end{bmatrix}$本质就是将三维基向量进行变形，使其x轴变换到$(\vec{x_x,x_y,x_z})$，y轴变换到$(\vec{y_x,y_y,y_z})$，z轴变换到$(\vec{z_x,z_y,z_z})$。所以，上述矩阵可以写成：

$$\mathbf{M}=\begin{bmatrix}r_x & u_x & v_x & 0\\r_y & u_y & v_y & 0\\r_z & u_z & v_z & 0\\0 & 0 & 0 & 1\end{bmatrix}^{-1}T^{-1}(t)$$

由于$r,u,v$都已经normalize了，不包含缩放，所以上述左侧的矩阵就是一个旋转矩阵，其逆矩阵等于转置矩阵，所以：

$$\mathbf{M}=\begin{bmatrix}r_x & r_y & r_z & 0\\u_x & u_y & u_z & 0\\v_x & v_y & v_z & 0\\0 & 0 & 0 & 1\end{bmatrix}\mathbf{T}(-\mathbf{t})=\begin{bmatrix}r_x & r_y & r_z & 0\\u_x & u_y & u_z & 0\\v_x & v_y & v_z & 0\\0 & 0 & 0 & 1\end{bmatrix}\begin{bmatrix}1 & 0 & 0 & -t_x\\0 & 1 & 0 & -t_y\\0 & 0 & 1 & -t_z\\0 & 0 & 0 & 1\end{bmatrix}\\=\begin{bmatrix}r_x & r_y & r_z & -(r_xt_x+r_yt_y+r_zt_z)\\u_x & u_y & u_z & -(u_xt_x+u_yt_y+u_zt_z)\\v_x & v_y & v_z & -(v_xt_x+v_yt_y+v_zt_z)\\0 & 0 & 0 & 1\end{bmatrix}=\begin{bmatrix}r_x & r_y & r_z & -\mathbf{t}\cdot \mathbf{r}\\u_x & u_y & u_z & -\mathbf{t}\cdot \mathbf{u}\\v_x & v_y & v_z & -\mathbf{t}\cdot \mathbf{v}\\0 & 0 & 0 & 1\end{bmatrix}$$

法线转换(Normal Transform)

本章中提到的转换矩阵都不适用于法线，下图就是将缩放矩阵应用于法线的错误例子：

正确的计算方式是乘以矩阵的伴随矩阵(adjoint)的转置矩阵，再将其normalize。

伴随矩阵：$\mathbf{A}$的伴随矩阵记作$adj(\mathbf{A})$，或$\mathbf{A}^*$，如果$\mathbf{A}$可逆，则$\mathbf{A}^*=\mathbf{A}^{-1}*det(\mathbf{A})$。

计算一个4x4矩阵的伴随矩阵的性能开销很大也没必要。位移并不会改变法线向量，而且大多数模型变换都是仿射的，并不会进行投影改变齐次坐标的w分量。因此在大多数情况下，都只需要计算左上角3x3矩阵的伴随矩阵即可。

通常甚至连伴随矩阵都不需要计算。假设一个变换矩阵只包含了位移、旋转和统一缩放。其中位移不会影响法线向量，统一缩放只会改变法线向量的长度，而旋转也只是改变了一个角度而已。传统的计算法线变换的方式就是求取变换矩阵的逆矩阵的转置矩阵，而旋转矩阵的逆矩阵和转置矩阵相等，它的逆矩阵的转置矩阵就是它自己。因此在这种情况下，模型的变换矩阵就能用来计算法线变换矩阵。

最后，normalize也并不是每次都需要执行的。假如变换只包含了位移和旋转，那法线向量的长度并没有被改变，因此不需要normalize。如果是统一缩放，那归一化时也只需要除以缩放系数就行，或者在计算时先将变换矩阵左上角3x3的部分除以缩放系数。

计算逆矩阵

下面三种方法可以计算逆矩阵：

• 如果由一个或一组简单变换组成的变换矩阵，则将其系数取负、顺序颠倒即可。如$\mathbf{M}=\mathbf{T}(\mathbf{t})\mathbf{R}(\phi)$，则$\mathbf{M}^{-1}=\mathbf{R}(-\phi)\mathbf{T}(\mathbf{-t})$
• 如果矩阵是正交矩阵，则逆矩阵为转置矩阵。一个或一组旋转矩阵都是正交矩阵。
• 如果不确定，则使用伴随矩阵、克莱默法则(Cramer’s Rule)、LU分解(LU Decomposition)、高斯消元法(Gaussian elimination)等方法计算逆矩阵。前两者由于使用更少的if分支语句，因此在现代GPU架构下更受欢迎。

克莱默法则：详细的推导和几何理解可以参考该视频。结论是：

如果$\mathbf{A}=\begin{bmatrix}a & b\\c & d\end{bmatrix},\mathbf{A}\cdot \begin{bmatrix}x\\y\end{bmatrix}=\begin{bmatrix}1\\2\end{bmatrix}$

则$x=\frac{det(\begin{bmatrix}1 & b\\ 2 & d\end{bmatrix})}{det(\mathbf{A})},y=\frac{det(\begin{bmatrix}a & 1\\ c & 2\end{bmatrix})}{det(\mathbf{A})}$

对于逆矩阵，我们知道$\mathbf{A}\mathbf{A}^{-1}=\mathbf{I}$。

假设$\mathbf{A}=\begin{bmatrix}a & b\\c & d\end{bmatrix},\mathbf{A}^{-1}=\begin{bmatrix}x_1 & x_2\\y_1 & y_2\end{bmatrix}$

则有：$\begin{bmatrix}a & b\\c & d\end{bmatrix}\begin{bmatrix}x_1\\x_2\end{bmatrix}=\begin{bmatrix}1\\0\end{bmatrix},\begin{bmatrix}a & b\\c & d\end{bmatrix}\begin{bmatrix}y_1\\y_2\end{bmatrix}=\begin{bmatrix}0\\1\end{bmatrix}$

根据克莱默法则可得： $x_1=\frac{d}{det(\mathbf{A})},x_2=\frac{-c}{det(\mathbf{A})},y_1=\frac{-b}{det(\mathbf{A})},y_2=\frac{a}{det(\mathbf{A})}$

因此：$\mathbf{A}^{-1}=\frac{1}{det(\mathbf{A})}\begin{bmatrix}d & -b\\-c & a\end{bmatrix}$

特殊矩阵变换和运算

欧拉变换

这个变换用来改变物体的朝向。

特殊矩阵转换和运算

欧拉转换(The Euler Transform)

用来转换物体的朝向，其矩阵可以有很多种顺序组合，但通常表示为：

$$\mathbf{E}(h,p,r)=\mathbf{R}_z(r)\mathbf{R}_x(p)\mathbf{R}_y(h)$$

由于E是多个旋转矩阵的组合，因此它是正交矩阵。它的逆矩阵$\mathbf{E}^{-1}=\mathbf{E}^{T}=(\mathbf{R}_z\mathbf{R}_x\mathbf{R}_y)^{T}=\mathbf{R}_y^{T}\mathbf{R}_x^{T}\mathbf{R}_z^{T}$。当然，更简单的方法是使用E自己的转置矩阵。

欧拉角(Euler angles)有时候也会被表示为rolls，即head=y-roll，pitch是x-roll。在飞行模拟中，head还会被表示为yaw。

不同的软件对于朝上的轴向定义不同，主要有z、-z、y轴朝向。本章节我们默认y轴朝向。

Euler angles适用于小角度转换或视角转向，但局限性也很大。比如，它很难处理两组E的组合，因为在两组E之间进行插值并不是简单的对h、p、r做插值，假设有两组完全不同的r、p、y，其实最后的结果E完全相同，那它们之间的插值也应该都是E，而如果对h、p、r做插值，得到的结果就不会和E相同。

从欧拉转换中提取参数(Extracting Parameters from the Euler Transform)

从E中推导出h、p、r三个值的过程如下（由于都是旋转矩阵，因此不需要使用4x4矩阵，只需要3x3矩阵）：

$$\mathbf{E}(h,p,r)=\begin{bmatrix}e_{00} & e_{01} & e_{02}\\e_{10} & e_{11} & e_{12}\\e_{20} & e_{21} & e_{22}\end{bmatrix}=\mathbf{R}_z(r)\mathbf{R}_x(p)\mathbf{R}_y(h)$$

$$=\begin{bmatrix}\cos r & -\sin r & 0\\\sin r & \cos r & 0\\0 & 0 & 1\end{bmatrix}\begin{bmatrix}1 & 0 & 0\\0 & cos p & -\sin p\\0 & \sin p & \cos p\end{bmatrix}\begin{bmatrix}\cos h & 0 & \sin h\\0 & 1 & 0\\-\sin h & 0 & \cos h\end{bmatrix}$$

$$=\begin{bmatrix}\cos r & -\sin r & 0\\\sin r & \cos r & 0\\0 & 0 & 1\end{bmatrix}\begin{bmatrix}\cos h & 0 & \sin h\\\sin p\sin h & \cos p & -\sin p\cos h\\-\cos p\sin h & \sin p & \cos p\cos h\end{bmatrix}$$

$$=\begin{bmatrix}\cos r\cos h-\sin r\sin p\sin h & -\sin r\cos p & \cos r\sin h+\sin r\sin p\cos h\\\sin r\cos h+\cos r\sin p\sin h & \cos r\cos p & \sin r\sin h-\cos r\sin p\cos h\\-\cos p\sin h & \sin p & \cos p\cos h\end{bmatrix}$$

$$\therefore e_{21}=\sin p,\;\;\;\frac{e_{01}}{e_{11}}=\frac{-\sin r}{\cos r}=-\tan r,\;\;\;\frac{e_{20}}{e_{22}}=\frac{-\sin h}{\cos h}=-\tan h$$

$$\therefore \begin{matrix}h=\arctan(-\frac{e_{20}}{e_{22}})=\text{atan2}(-e_{20},e_{22})\\p=\arcsin(e_{21})\\r=\arctan(-\frac{e_{01}}{e_{11}})=\text{atan2}(-e_{01},e_{11})\end{matrix}$$

上述推导还可能遇到$\cos p=0$的特殊情况。此时r、h都会使物体沿着同一个轴旋转，这种情况叫做万向锁(gimbal lock)(更详细说明参考该视频)。在这种情况下，r的值决定了在这个轴上的最终旋转。因此，我们将h设为0，则：

$$\mathbf{E}=\begin{bmatrix}\cos r & \sin r\cos p & \sin r\sin p\\\sin r & \cos r\cos p & -\cos r\sin p\\0 & \sin p & \cos p\end{bmatrix}$$

$$\therefore \frac{e_{10}}{e_{00}}=\frac{\sin r}{\cos r}=\tan r$$

$$\therefore r=\text{atan2}(e_{10},e_{00})$$

需要注意的是，由于$p=\text{arcsin}(e_{21})$，因此$-\pi/2\leq p \leq \pi/2$。显然，如果构建E的时候p的值不在这个范围内，那根据上述步骤得到的p肯定和实际的p不一样。这也就意味着不止一组r、p、h可以得到同样的欧拉变换。

尽管欧拉角存在万向锁这种无法避免的问题，它依然被广泛使用，因为它在使用动画曲线编辑器时很好用。

矩阵分解(Matrix Decomposition)

前面的例子都是假设我们已经知道物体经过了哪些变换得到了最终的变换矩阵。但实际运用时我们往往不知道中间经过了哪些步骤，这就需要分解矩阵。

平移矩阵可以直接通过变换矩阵的最后一列获取；计算行列式的正负可以判断是否发生了反射；旋转、缩放和剪切不是这么直观就能分解出来的，但是也有一些固定的方法，但本文不再讲述。

绕任意轴旋转(Rotate about an Arbitrary Axis)

假设我们需要将物体绕r轴旋转ɑ度，我们需要先对物体进行M变换，使其本地空间的x轴和r轴对齐，将其旋转后再进行M^-1变换，如下图所示。

通过观察我们可以得知：

$$\mathbf{s}=\begin{cases} (0,-r_z,r_y) &\text{ if } |r_x| \leq |r_y| \;\;\;\text{and}\;\;\; |r_x| \leq |r_z| \\ (-r_z,0,r_x) &\text{ if } |r_y| \leq |r_x| \;\;\;\text{and}\;\;\; |r_y| \leq |r_z| \\ (-r_y,r_x,0) &\text{ if } |r_z| \leq |r_x| \;\;\;\text{and}\;\;\; |r_z| \leq |r_y| \end{cases}$$

将s和r normalize之后，$\mathbf{t}=\mathbf{s} \times \mathbf{r}$。还有一些其他方法可以省略计算s的分支步骤。将r、s、t和x、y、z轴对齐的转换矩阵为：

$$\mathbf{M}=\begin{bmatrix}\mathbf{r}^T\\\mathbf{s}^T\\\mathbf{t}^T\end{bmatrix}$$

最终用来进行转换的矩阵为：

$$\mathbf{X}=\mathbf{M}^T\mathbf{R}_x(\alpha)\mathbf{M}$$

还有另外一个直接进行旋转变换的公式将不做推导，直接给出，如下：

$$\tiny{ \mathbf{R}=\\ \begin{bmatrix} \cos \phi+(1-\cos \phi)r_x^2 & (1-\cos \phi)r_xr_y-r_z\sin \phi & (1-\cos \phi)r_xr_z+r_y\sin \phi\\ (1-\cos \phi)r_xr_y+r_z\sin \phi & \cos \phi+(1-\cos \phi)r_y^2 & (1-\cos \phi)r_yr_z-r_x\sin \phi\\ (1-\cos \phi)r_xr_z-r_y\sin \phi & (1-\cos \phi)r_yr_z+r_x\sin \phi & \cos \phi+(1-\cos \phi)r_z^2 \end{bmatrix} }$$

四元数(Quaternions)

虚数(imaginary number)：通俗表示为$i=\sqrt{-1}$。实数乘以i可以理解为将一维的实数域拓展到了二维的复数域，将一维上的实数逆时针旋转90°，$i^2$则是旋转180°，也即$1\times i \times i=-1$，即正数变成负数。 复数(complex number)：表示为$Z=a+bi$，a、b为实数，i为虚数，a称为实部，b称为虚部。在复数域上，它代表一个位于(a,b)的点(数字)。

数学背景(Mathematical Background)

我们将四元数$\hat{\mathbf{q}}$定义为：

$$\hat{\mathbf{q}}=(\mathbf{q}_v,q_w)=iq_x+jq_y+kq_z+q_w=\mathbf{q}_v+q_w$$

$$\mathbf{q}_v=iq_x+jq_y+kq_z=(q_x,q_y,q_z)$$

$$i^2=j^2=k^2=-1,jk=-kj=i,ki=-ik=j,ji=-ji=k$$

其中，$q_w$被称为$\hat{\mathbf{q}}$的实部，$\mathbf{q}_v$被称为虚部，i、j、k为虚数单位。$\mathbf{q}_v$适用所有法向量运算，如加法、缩放、点乘、叉乘等。两个四元数$\hat{\mathbf{q}},\hat{\mathbf{r}}$相乘(顺序不可逆)的推导过程如下：

$$\hat{\mathbf{q}}\hat{\mathbf{r}}=(iq_x+jq_y+kq_z+q_w)(ir_x+jr_y+kr_z+r_w)$$

$$=i(q_yr_z-q_zr_y+r_wq_x+q_wr_x)$$

$$+j(q_zr_x-q_xr_z+r_wq_y+q_wr_y)$$

$$+k(q_xr_y-q_yr_x+r_wq_z+q_wr_z)$$

$$+q_wr_w-q_xr_x-q_yr_y-q_zr_z$$

$$=(\mathbf{q}_v\times \mathbf{r}_v+r_w\mathbf{q}_v+q_w\mathbf{r}_v,q_wr_w-\mathbf{q}_v\cdot\mathbf{r}_v)$$

其加法(addition)、共轭(conjugate)、模(norm)、恒等式(Identity)分别为：

$$\textbf{Addition:}\;\;\;\hat{\mathbf{q}}+\hat{\mathbf{r}}=(\mathbf{q}_v,q_w)+(\mathbf{r}_v,r_w)=(\mathbf{q}_v+\mathbf{r}_v,q_w+r_w)$$

$$\textbf{Conjugate:}\;\;\;\hat{\mathbf{q}}^*=(\mathbf{q}_v,q_w)^*=(-\mathbf{q}_v,q_w)$$

$$\textbf{Norm:}\;\;\;n(\hat{\mathbf{q}})=\|\hat{\mathbf{q}}\|=\sqrt{\hat{\mathbf{q}}\hat{\mathbf{q}}^*}=\sqrt{\hat{\mathbf{q}}^*\hat{\mathbf{q}}}=\sqrt{\mathbf{q}_v\cdot \mathbf{q}_v+q_w^2}$$

$$=\sqrt{q_x^2+q_y^2+q_z^2+q_w^2}$$

$$\textbf{Identity:}\;\;\;\hat{\mathbf{i}}=(\mathbf{0},1)$$

由于$\hat{\mathbf{q}}^{-1}\hat{\mathbf{q}}=\hat{\mathbf{q}}\hat{\mathbf{q}}^{-1}=1$，因此：

$$\hat{\mathbf{q}}^{-1}=\frac{1}{\hat{\mathbf{q}}}=\frac{1}{\|\hat{\mathbf{q}}\|^2}\hat{\mathbf{q}}^*$$

四元数乘以一个实数的运算和顺序无关，即：

$$s\hat{\mathbf{q}}=\hat{\mathbf{q}}s=(s\mathbf{q}_v,sq_w)$$

单位四元数(unit quaternion)即$\|\hat{\mathbf{q}}\|=1$。我们知道复数(单位向量)可以表示为：

$$\mathbf{Z}=i\sin \phi+\cos \phi=e^{i\phi}$$

我们知道$\frac{d}{d_x}e^{ix}=ie^{ix}$，即$f(x)=e^x$的导数是它自己，$f(x)=e^{ix}$的导数是i乘以它自己。对于复数$Z=f(x)=i\sin x+\cos x$，它的导数其实就是这个点在复平面上沿着单位圆的切线，也就是这个点乘以虚数i。因此，$\frac{d}{d_x}f(x)=if(x)$，根据e的定义，可以判定$f(x)=e^{ix}$

类似的，四元数也能表示为：

$$\hat{\mathbf{q}}=(\sin \phi \mathbf{u}_q,\cos \phi)=\sin \phi \mathbf{u}_q+\cos \phi$$

其中，$\mathbf{u}_q$为三维向量，且$\|\mathbf{u}_q\|=1$

四元数转换(Quaternion Transform)

我们可以很方便地用单位四元数来表示物体在三维空间的旋转。假设现在有一个点或向量$\mathbf{p}=\begin{bmatrix}p_x & p_y & p_z & p_w\end{bmatrix}^T$，我们用它的四个坐标分量作为四元数$\hat{\mathbf{p}}$的四个实数参数。另有一个单位四元数$\hat{\mathbf{q}}=(\sin \phi \mathbf{u}_q,\cos \phi)$。可以证明：$\hat{\mathbf{q}}\hat{\mathbf{p}}\hat{\mathbf{q}}^{-1}$将$\hat{\mathbf{p}}$和p沿着$\mathbf{u}_q$旋转了$2\phi$角度。由于$\hat{\mathbf{q}}$是单位四元数，因此$\|q\|=1$。由4.3.1的公式可以得出，$\hat{\mathbf{q}}^{-1}=\hat{\mathbf{q}}^*$。

给定两个单位四元数$\hat{\mathbf{q}},\hat{\mathbf{r}}$，对p先按$\hat{\mathbf{q}}$变换再按$\hat{\mathbf{r}}$变换的公式为：

$$\hat{\mathbf{r}}(\hat{\mathbf{q}}\hat{\mathbf{p}}\hat{\mathbf{q}}^*)\hat{\mathbf{r}}^*=(\hat{\mathbf{r}}\hat{\mathbf{q}})\hat{\mathbf{p}}(\hat{\mathbf{r}}\hat{\mathbf{q}})^*=\hat{\mathbf{c}}\hat{\mathbf{p}}\hat{\mathbf{c}}^*$$

其中，$\hat{\mathbf{c}}=\hat{\mathbf{r}}\hat{\mathbf{q}}$。

本段后续内容都是对于各种变换的计算公式，由于这些在软件或游戏引擎的底层中都已经实现，我们大概了解原理即可，具体推导过程省略。

顶点混合(Vertex Blending)

假如现在有一支机械臂，它的动画由上臂和前臂组成。如果这个动画通过刚体变换来实现，就显得不是那么真实，因为关节包含了上臂和前臂重叠的部分。更好的实现方式是将整个机械臂当成一个物体来处理。

Vertex Blending又称为线性混合蒙皮(linear-blend skinning)、包络(enveloping)或骨骼子空间变形(skeleton-subspace deformation)。简单来说，在上述例子中，前臂和上臂各自使用单独的动画，但是在关节处通过一个拥有弹性的部位相连。这个部位有一部分顶点受上臂的变换矩阵影响，有一部分顶点受前臂的变换矩阵影响。在此基础上，同一个顶点就可以受到多个变换矩阵的影响，最终将不同权重的变换矩阵结果混合起来。这种用户自定义某个顶点受哪些骨骼影响的权重的操作就是所谓的“刷权重”。一个顶点p随着时间t变换的坐标可以表示为：

$$\mathbf{u}(t)=\sum_{i=0}^{n-1}w_i\mathbf{B}_i(t)\mathbf{M}_i^{-1}\mathbf{p}$$

其中，n为所有骨骼的数量，$w_i$为某个骨骼对于顶点p的权重。$\mathbf{M}_i$将骨骼从本地坐标系转换到世界坐标系，$\mathbf{B}_i(t)$为骨骼在当前时间点在世界坐标系中的转换矩阵。

Vertex Blending很适合在GPU上使用。网格中的顶点可以存储在GPU的静态缓冲区，每一帧只需要传入骨骼的变换矩阵，就能通过顶点着色器计算出所有顶点所受的影响，这样可以有效降低CPU和GPU之间的数据传输，提高GPU渲染效率。骨骼变换的信息也可以存储在贴图中，以防寄存器存储达到上限。

在一些特殊的算法，比如下一节的变形目标(morph targets)中，所有的权重加起来也可以不为1。

线性的Vertex Blending算法的缺点是容易产生不必要的折叠(folding)、扭曲(twisting)和子相交(self-intersection)，如下图所示。

更好的办法是使用对偶四元数(dual quaternions)。这种方法可以在蒙皮时保持原始变形的刚性，避免产生想糖纸一样的扭曲。和linear blend skinning相比，这种算法消耗仅为1.5倍不到，因此很快就流行起来。但是dual quaternions会导致肿胀(bulging)效应，因此又出现了旋转中心蒙皮(center-of-rotation skinning)算法。这种算法假设：

1. 局部变换应该是刚体变换
2. 拥有相似权重的顶点应该具有相似的变换

每个顶点的旋转中心都是预先计算好的，正交(刚体)约束则用来防止弯头坍塌和糖纸扭曲。在运行时，该算法类似于linear blend skinning，旋转中心由GPU执行linear blend skinning，然后执行四元数混合步骤。

目前最常用的算法似乎还是Linear Blend Skinning(LBS)。Center-of-Rotation Skinning(CoRS)出自Disney的一篇论文，从其展示的效果和效率来看，似乎是目前效果最好、效率相对较高的算法。

变形(Morphing)

简单来说就是在时间点t0和t1分别有两个不同的三维模型，通过某种插值算法在两个时间点内的每个时间点都能获得一个模型。

Morphing主要需要解决两个问题，一个是顶点对应(vertex correspondence)，一个是插值(interpolation)。给定两个任意模型，它们可能拥有不同的拓扑结构、顶点数量和网格连接方式，因此通常需要先预设它们的顶点对应关系。这是一个难题，已经有很多相关研究。

如果两个模型已经拥有一对一的顶点关系，那就可以很容易地对每个顶点进行插值（比如线性插值）。

有一种变形方式叫做变形目标(morph targets)或混合形状(blend shapes)，用户可以对变形进行更直观的控制。

在上图中，我们有一张中性的脸的模型，记作$\mathcal{N}$。它可以有一系列($k\geq 1$个)表情，比如图中的笑脸。我们将这一系列表情记作$\mathcal{P}_i,i\in [1,…,k]$，它们和中性脸的差异记作$\mathcal{D}_i=\mathcal{P}_i-\mathcal{N}$。因此，一张变形的脸的模型可以表示为：

$$\mathcal{M}=\mathcal{N}+\sum_{i=1}^{k}w_i\mathcal{D}_i$$

其中$w_i$为每种差异的权重，如果为1，则得到了上图正中间的笑脸，如果为0.5，则得到了一张半笑的脸。当然，权重也可以小于0(伤心的脸)或大于1(狂笑的脸)。

Morph targets是动画中很有用的一项技术，因为用户可以独立控制每个不同的动画特征。pose-space deformation技术可以将vertex blending和morph targets结合起来。还可以使用贴图来存储和读取动画的切换。支持stream-out和每个顶点ID的硬件允许单个模型使用更多数量的变形目标，并且可以专门在GPU上进行计算。使用低分辨率的网格、在曲面细分阶段生成更多顶点并使用位移贴图(displacement texture)可以避免给高分辨率网格的每个顶点进行蒙皮带来的消耗。

几何缓存回放(Geometry Cache Playback)

在剪辑场景中，可能需要使用高质量的动画，这种变形用上面的算法无法表示出来。一个很自然的做法是将每一帧所有的顶点数据存储在本地磁盘，实时读取这些数据来更新网格。但是对于一段短动画中一个简单的30000个顶点的模型，这些数据的传输量就能达到50MB/s。有一些方法可以将内存消耗降低到10%左右。

首先，使用量化(quantization)。比如，对每个坐标用16位整数存储其位置和贴图坐标。这个步骤是有损的，在压缩数据后无法还原原始数据。为了进一步减少数据量，会进行空间(spatial)和时间(temporal)上的预测并记录插值。对于空间压缩，可以使用平行四边形预测(parallelogram prediction)：对于一个三角形带，下一个顶点的预测位置应该与当前三角形形成一个平行四边形，记录下这个预测位置和实际位置的差。如果预测理想的话，这个差值应该接近于0，大多数压缩算法都能对其进行很好地压缩。对于时间压缩，做法和MPEG类似：每隔n帧进行时间压缩。如果一个顶点在第n-1帧和第n帧之间移动了一个向量，那么它在第n帧和第n+1帧之间也应该移动类似的向量。这些技术能显著降低数据存储量以进行实时串流。

投影(Projection)

在本章前面提到的变换都没有涉及第四个坐标w分量，在变换后，点还是点，向量还是向量，4x4矩阵的最后一行依然是[0 0 0 1]。透视投影矩阵(Perspective Projection Matrices)却不一样，它的w分量往往不是1，其他分量需要除以w。4.7.1中介绍的正交投影(Orthographic Projection)则和w分量无关。

在本节中，我们默认观察者的视线是沿着z轴的负方向，上方为y轴，右手为x轴。

正交投影(Orthographic Projection)

正交投影的一个特征是平行线投影完依然是平行线。不管物体距离相机多远，它的大小都不变，它的变换矩阵可以表示为：

$$\mathbf{P}_o=\begin{bmatrix} 1 & 0 & 0 & 0\\ 0 & 1 & 0 & 0\\ 0 & 0 & 0 & 0\\ 0 & 0 & 0 & 1 \end{bmatrix}$$

由于$\|\mathbf{P}_o\|=0$，因此$\mathbf{P}_o$不可逆，换句话说，从三维投影到二维后，无法恢复到三维。由于-z轴上的点和z轴上的点都会被投影到二维平面上，我们通常都会将z轴的投影范围限制在一定范围内，即从n(near plane/front plane/hither)到f(far plane/back plane/yon)。

另一种更常用的正交投影矩阵可以表示为$[l, r, b, t, n, f]$，即左、右、底、顶、近、远平面。这个矩阵将由这6个平面围成的轴对齐包围盒(axis-aligned bounding box, AABB)变换成一个以原点为中心的轴对齐立方体，如下图所示：

AABB的最小坐标为(l, b, n)，最大坐标为(r, t, f)。由于观察者是看向z轴负方向，因此n>f。在OpenGL中，变换后的轴对齐立方体的最小坐标为(-1, -1, -1)，最大坐标为(1, 1, 1)。在DirectX中两者则是(-1, -1, 0)和(1, 1, 1)。这个立方体就是第二章中提到的规范视域体(canonical view volume)，与之对应的坐标就是标准化设备坐标(Normalized Device Coordinate, NDC)。这个变换的作用是方便后续的裁剪阶段。

完整的投影变换可以表示为：

$$\mathbf{P}_o=\mathbf{S}(\mathbf{s})\mathbf{T}(\mathbf{t})\\ =\begin{bmatrix} \frac{2}{r-l} & 0 & 0 & 0\\ 0 & \frac{2}{t-b} & 0 & 0\\ 0 & 0 & \frac{2}{f-n} & 0\\ 0 & 0 & 0 & 1 \end{bmatrix}\begin{bmatrix} 1 & 0 & 0 & -\frac{r+l}{2}\\ 0 & 1 & 0 & -\frac{t+b}{2}\\ 0 & 0 & 1 & -\frac{f+n}{2}\\ 0 & 0 & 0 & 1 \end{bmatrix}\\ =\begin{bmatrix} \frac{2}{r-l} & 0 & 0 & -\frac{r+l}{r-l}\\ 0 & \frac{2}{t-b} & 0 & -\frac{t+b}{t-b}\\ 0 & 0 & \frac{2}{f-n} & -\frac{f+n}{f-n}\\ 0 & 0 & 0 & 1 \end{bmatrix}$$

其中，$\mathbf{s}=(2/(r-l),2/(t-b),2/(f-n)),\mathbf{t}=(-(r+l)/2,-(t+b)/2,-(f+n)/2)$。$\mathbf{P}_o$是可逆的，$\mathbf{P}_o^{-1}=\mathbf{T}(-t)\mathbf{S}((r-l)/2,(t-b)/2,(f-n)/2)$