这两节内容是关于光栅化，即在我们通过前两节内容得到观测变换后的模型信息后，通过光栅化的方式把图像画在屏幕上。

Lecture 5 光栅化【三角形的离散化】

在 MVP 变换后，我们将模型缩到了 $[-1,1]^3$ 大小的空间中，光栅化的目标就是将他们显示在屏幕上。我准备以几个步骤来构建整个过程。

光栅化 Rasterize 是来自德语，Raster 是德语的屏幕，Rasterize 表示在屏幕上绘制的意思

一、构建一个屏幕

我们都知道屏幕是一个平面，上面有多个像素点，那我们可以把每个像素值都当作是一个二维数组中的值，用一个分辨率大小的二维数组来描述整个屏幕，就像下面这个图：

那么很自然能想到，我们在访问的时候就是这样的操作：

for (int x = 0; x < width; x++)
{
  for (int y = 0; y < height; y++)
  {
    // 俺就是像素值
    cout << screen[x][y] << endl;
  }
}

现在我们可以访问屏幕上的像素值了，但是我们的模型所在的空间是 $[-1,1]^3$，屏幕却是 $[0, width] \times [0, height]$，因此需要一个视口变换来将模型转换到屏幕空间上。

二、视口变换

屏幕是没有 $z$ 轴的，因此我们暂时不管 $z$ 轴上的点，不过后面在判断物体在空间中深度位置的时候还是会用到的。

• 视口变换矩阵：

  $$M_{viewport} = \begin{pmatrix} \frac{width}{2} & 0 & 0 & \frac{width}{2} \\ 0 & \frac{height}{2} & 0 & \frac{height}{2} \\ 0 & 0 & 1 & 0 \\ 0 & 0 & 0 & 1 \end{pmatrix}$$

三、三角形离散化表示

我们通过变换得到屏幕空间中的模型后，需要在渲染的时候定义一个图元结构（Primitive），用图元结构来离散化组成整个图像，一般有：三角形、四边形、多边形、线段、点

为什么选择三角形呢？

• 三角形是最基础的几何面，可以组成任意多边形（你应该是找不到比这个更简单的了）
• 在三维空间中，三角形 3 个顶点组成的面肯定是一个平面，4 个顶点组成的可不一定是一个四边形，也可能直接变成一个锥体
• 三角形一定是凸多边形（凹多边形容易导致算法错误，或者在算法实现上更困难）
• 三角形可以用更方便的方法做插值（比如利用重心坐标做插值）

还有什么理由不选它呢，是吧。

然后我们选择了三角形来作为图像的表示，但是在屏幕上显示的时候，像素是有大小的一种离散化表示，于是我们得把三角形也离散化，最简单的离散化近似方式就是采样。

就像上图，我们取每一个像素的中心点，如果中心点在三角形内，那就“点亮”这个像素，如果不在，就不亮。因此光栅化过程可以写成这样：

for (int x = 0; x < width; x++)
{
  for (int y = 0; y < height; y++)
  {
    // 加 0.5 是取中心的坐标
    image[x][y] = inside(tri, x + 0.5, y + 0.5);
  }
}

这里还是有一个问题，怎么求点在三角形中呢，最简单的方法就是用叉乘的性质：

1. 我们知道叉乘的结果是一个垂直于两个向量的向量
2. 用顶点 A 指向顶点 B 的向量 $\vec{AB}$ 叉乘顶点 A 到像素中点 P 的向量 $\vec{AP}$ 可以得到另一个向量 $\vec{c}$
3. 通过右手定则可以由向量 $\vec{c}$ 的方向判断出点 p 在 $\vec{AB}$ 的左边还是右边
4. 对三角形的每个边按特定模式都进行步骤 3，就可以判断出到底点在三角形内还是外

就比如这个图，我们求得 $\vec{P_2P_0}$ 和 $\vec{P_2Q}$，通过 $\vec{P_2P_0} \times \vec{P_2Q}$ 可得到 $\vec{c_1}$ 朝向屏幕内，通过这样的模式，可得到 Q 在 $\vec{P_2P_0}$ 右边，在 $\vec{P_0P_1}$ 左边，在 $\vec{P_1P_2}$ 左边，于是就可以知道 Q 点在三角形外。

四、优化速度

基于包围盒的采样

对每个像素都进行采样的话是非常浪费资源的，比如上图，很明显最左边的那一列是完全不需要纳入判断的，只有 $[x_{p_0}, x_{p_1}] \times [y_{p_0}, y_{p_2}]$ 范围内的点才有可能落在三角形内，所以只需要遍历这个区域，我们称这样四个顶点的最大最小值围成的矩形为 Bounding Box。

这是一个最优解了吗，其实并不是，当三角形过于细长但是角度又比较偏的时候，三角形的面积占比只会占包围盒的一小部分。

基于自增的采样

对于倾斜三角形，我们可以从每行的一边开始从左到右进行采样：

1. 从顶点 x 和 y 值最小的一个点开始作为起点
2. 从左往右进行遍历，遇到点在三角形外的情况就停止
3. 回到起点，然后往上一行，重新开始寻找新起点，遇到点在三角形内的情况停止
4. 找到新起点后重复 2，3 步骤，直到找不到新起点（最后一行会遍历到最右边）

这种方法的循环次数会少一些，但是实现起来比较复杂，而且不容易进行并行计算（显卡里大部分这种过程是并行计算，所以会很快），因此实际情况可能还会慢。

五、锯齿问题

锯齿问题在光栅化过程中是不可避免的，但是我们可以通过一些反走样的方法让锯齿看起来并不是那么严重。

Lecture 6 光栅化【深度测试与抗锯齿】

这里我们先明确一个概念，光栅化是对 2D 点的采样，照片是对相机传感器平面进行采样，视频是对时间进行采样（帧）。

反走样【Anti-Aliasing】

走样出现的原因：信号变化太快【高频】，而采样的速度太慢

常见锯齿类型：

• 光栅化：

  

• 拍照：

  

• 视频：

  

详细解释走样问题之前，我们需要先回到信号学的角度才能理解其中的本质：信号本身就是一个函数，输入参数，输出结果。对于任何的函数 $f(x)$，都有两种表现形式：时域、频域。其中我们最常见的便是时域的形式，也是现实世界的直观形式，可以理解成即 $f(x)$ 本身，把 $x$ 当成时间，那么时域便是分析函数值关于时间 $x$ 的关系。

傅里叶告诉我们，任意一个函数 $f(x)$ 可以由若干频率不同，振幅不同的余弦函数相加而成。
详见：为什么傅里叶变换可以把时域信号变为频域信号？ - 李狗嗨的回答 - 知乎

当我们每隔一段时间采样 $f(x)$，就可以看作是每隔一段时间，采样组成 $f(x)$ 的一系列余弦波然后相加起来。如下图所示，低频采样的结果连接起来的曲线和真实值是比较接近的，但是在高频区域，信号的采集连接起来的曲线和真实情况就出现了巨大的失真。

就像上图表现的一样，当我们用白点的频率采样蓝色曲线的时候，最后采样出来的结果拟合出来的黑色曲线看起来就像是一个低频信息，根本无法用这个结果对真实信息进行区分，所以也叫走样。

解决的方法：过滤掉部分频率的内容

当我们通过傅里叶变换将一幅图像从时域转换到频域上，就能得到右图，这个时候我们会发现，大多数图像的信息是处于中心位置的，也就是低频区域。

• 对图像的低频信息进行过滤（高通滤波）

• 对图像的高频信息进行过滤（低通滤波）

在卷积理论中我们知道，在时域上对图像进行卷积的结果，等同于在频域上图像与滤波的乘积做逆傅里叶变换。（时域的卷积 = 频域的乘积 || 频域的卷积 = 时域的乘积）

从频域上理解采样：采样就是重复频率上的内容

这段理解起来有点复杂，我尽量记一下，如果看不懂可以去听课里这段。

从上图中可以看出，采样相当于是用狄卡拉梳状函数（冲激函数）乘以时域上的函数，相当于是在频域上卷积，最后得到重复的一堆函数。知乎上有一篇文章对这个过程进行了详细解释。

因此我们可以把采样的走样情况理解为频域上内容的混叠，如下图：

在 $\pm\frac{F_{s_d}}{2}$ 处，当我们进行密集采样的时候，生成的频谱图像可以看出来，每个信号互相并不会有干扰发生，也就是没有混叠的情况，但是当我们降低采样率的时候，频谱图像就发生了混叠，也就产生了走样现象。因此也可以从这个角度来解释反走样问题。

对图像进行模糊操作相当于过滤掉了高频部分，对应到频谱图像上就是将信号在宽度上缩减避免了混叠的情况出现。

上面的操作转换时域上就是上图这样，看起来效果还不错。

除了这种简单的反走样方法，我们现在普遍使用的是通过超采样的方法来让图像变得模糊。即下面的 MSAA（Multisample anti aliasing）：

原本图像是对单个像素点判断是否在三角形内还是外来实现光栅化过程，但是在 MSAA 中，我们将单个像素细分为 $2 \times 2$ 的“子像素”，通过对“子像素”在三角形中的比例进行像素值的模糊。

MSAA 并不是通过提升图像分辨率解决的反走样问题，只是一个将图像模糊的方式。（并没有提高对信号的采样率）

最后得到上图的结果，就是等于对一个像素内的多个位置进行采样并取其平均值来近似实现一个低通滤波器的效果。

MSAA 对计算量的消耗还是挺高的，在工业中，我们可以重新定义小像素的排列位置实现小像素在像素间的复用，因此实际上玩游戏的时候体验上并不会降低很多。

其他目前常用的抗锯齿方案

• FXAA(Fast Approximate AA)
• TAA(Temporal AA)
• DLSS(Deep Learning Super Sampling)

作业 2 题解

在作业 1 中，我们实现了在屏幕上画一个线框三角形，并且可以通过 A、D 键进行旋转，但是在作业 2 中我们要更进一步，在屏幕上画出一个实心三角形，也就是对一个三角形进行光栅化。

基础题

• 正确实现三角形栅格化算法
• 正确测试点是否在三角形内
• 正确实现 z-buffer 算法，将三角形按顺序画在屏幕上

提高题

• 用 Super-sampling 处理 Anti-aliasing（MSAAx2）

实现

想了一下，反正代码网上都挺多了，直接做一个代码记录吧，也比较省事。这次作业的目标是完成 rasterize_triangle(const Triangle& t) 函数并调用：

1. 创建三角形的 2 维 bounding box:

   // 通过找到在左下角一点的坐标和右上角一点的坐标来确定一个 bounding box
    float l_x = std::min(std::min(v[0].x(), v[1].x()), v[2].x());
    float l_y = std::min(std::min(v[0].y(), v[1].y()), v[2].y());
    float r_x = std::max(std::max(v[0].x(), v[1].x()), v[2].x());
    float r_y = std::max(std::max(v[0].y(), v[1].y()), v[2].y());
   // 四舍五入到整数范围进行索引
    l_x = (int)std::floor(l_x);
    l_y = (int)std::floor(l_y);
    r_x = (int)std::ceil(r_x);
    r_y = (int)std::ceil(r_y);

2. 遍历此 bounding box 内的所有像素（整数索引）。然后，使用像素中心的屏幕空间坐标来检查中心点是否在三角形内:

   // 判断是否在三角形内
    static bool insideTriangle(int x, int y, const Vector3f* _v)
    {
      const Eigen::Vector2f P(x, y);
      const Eigen::Vector2f A = _v[0].head(2), B = _v[1].head(2), C = _v[2].head(2);
   
      const Eigen::Vector2f AP = P - A;
      const Eigen::Vector2f AB = B - A;
      const Eigen::Vector2f BP = P - B;
      const Eigen::Vector2f BC = C - B;
      const Eigen::Vector2f CP = P - C;
      const Eigen::Vector2f CA = A - C;
   // 叉乘得到法向量，通过正负加上右手螺旋定则判断点在三角形内外
      float c1 = AB[0] * AP[1] - AB[1] * AP[0];
      float c2 = BC[0] * BP[1] - BC[1] * BP[0]; 
      float c3 = CA[0] * CP[1] - CA[1] * CP[0];
   
      if (c1 > 0 && c2 > 0 && c3 > 0) {
          return true;
      } else if (c1 < 0 && c2 < 0 && c3 < 0) { 
          return true;
      } else {
          return false;
      }
    }

3. 如果在内部，则将其位置处的插值深度值（interpolated depth value）与深度缓冲区（depth buffer）中的相应值进行比较，如果当前点更靠近相机，请设置像素颜色并更新深度缓冲区:

   for(int x = l_x; x < r_x; x++) {
     for(int y = l_y; y < r_y; y++) {
       if (insideTriangle(x + 0.5, y + 0.5, t.v)) {
         // computeBarycentric2D 是计算重心坐标，用于对后面对 z 值进行插值
         auto[alpha, beta, gamma] = computeBarycentric2D(x, y, t.v);
         float w_reciprocal = 1.0/(alpha / v[0].w() + beta / v[1].w() + gamma / v[2].w());
         float z_interpolated = alpha * v[0].z() / v[0].w() + beta * v[1].z() / v[1].w() + gamma * v[2].z() / v[2].w();
         z_interpolated *= w_reciprocal;
         // 深度测试
         if (z_interpolated < depth_buf[get_index(x, y)]) {
             set_pixel(Eigen::Vector3f(x, y, z_interpolated), t.getColor());
             depth_buf[get_index(x, y)] = z_interpolated;
         }
       }
     }
   }

   也许你也像我一样，看了很久 w() 这个是什么，这个其实是 eigen 中关于获取四元数中的 w 分量的函数，具体解释可以看这里：Eigen::Quaternion< Scalar, Options > Class Template Reference

4. 提高部分：

   bool MSAA = false;
   if (MSAA) {
       // 创建 2x2 的判断点
       std::vector<Eigen::Vector2f> pos
       {
           {0.25, 0.25},
           {0.75, 0.25},
           {0.25, 0.75},
           {0.75, 0.75}
       };
       for(int x = l_x; x < r_x; x++) {
           for(int y = l_y; y < r_y; y++) {
               float minDepth = FLT_MAX;
               int count = 0;
               for(int i = 0; i < 4; i++) {
                 // 其实也就是在这一步，将判断三角形内外的点设定为了 4 个点，也就完成了 MSAAx2
                   if (insideTriangle(x + pos[i][0], y + pos[i][1], t.v)) {
                       auto[alpha, beta, gamma] = computeBarycentric2D(x, y, t.v);
                       float w_reciprocal = 1.0/(alpha / v[0].w() + beta / v[1].w() + gamma / v[2].w());
                       float z_interpolated = alpha * v[0].z() / v[0].w() + beta * v[1].z() / v[1].w() + gamma * v[2].z() / v[2].w();
                       z_interpolated *= w_reciprocal;
                       minDepth = std::min(minDepth, z_interpolated);
                       count++;
                   }
               }
               if (count != 0) {
                   if (depth_buf[get_index(x, y)] > minDepth) {
                       Vector3f color = t.getColor() * count / 4.0;
                       Vector3f point(3);
                       point << (float)x, (float)y, minDepth;
                       // 替换深度
                       depth_buf[get_index(x, y)] = minDepth;
                       // 修改颜色
                       set_pixel(point, color);
                   }
               }
           }
       }
   } else {
       for(int x = l_x; x < r_x; x++) {
           for(int y = l_y; y < r_y; y++) {
               if (insideTriangle(x + 0.5, y + 0.5, t.v)) {
                   auto[alpha, beta, gamma] = computeBarycentric2D(x, y, t.v);
                   float w_reciprocal = 1.0/(alpha / v[0].w() + beta / v[1].w() + gamma / v[2].w());
                   float z_interpolated = alpha * v[0].z() / v[0].w() + beta * v[1].z() / v[1].w() + gamma * v[2].z() / v[2].w();
                   z_interpolated *= w_reciprocal;
                   if (z_interpolated < depth_buf[get_index(x, y)]) {
                       set_pixel(Eigen::Vector3f(x, y, z_interpolated), t.getColor());
                       depth_buf[get_index(x, y)] = z_interpolated;
                   }
               }
           }
       }
   }

   来看看具体效果吧，其实这次作业困扰了两天都无法理解是为什么，我先写了一个版本，然后发现没有效果，然后换成了别人的代码，也没有效果，对比了所有代码都是一模一样的，编译出来也没有效果，就比较离谱，然后直接用别人的代码编译就正常了，怪，太怪了。（Cmake 配置也是一样的都不行）

   

   

参考文章

[1] [实时渲染基础（2）光栅化（Rasterization）](https://www.cnblogs.com/KillerAery/p/14518859.html)
[2] [GAMES-101 现代计算机图形学入门笔记](https://blog.liukairui.me/article/GAMES-101%E7%8E%B0%E4%BB%A3%E8%AE%A1%E7%AE%97%E6%9C%BA%E5%9B%BE%E5%BD%A2%E5%AD%A6%E5%85%A5%E9%97%A8%E7%AC%94%E8%AE%B0/)
[3] [GAMES101-现代计算机图形学学习笔记（05）](https://blog.csdn.net/qq_36242312/article/details/105695242)
[4] [采样定理，频谱混叠和傅里叶变换 深入理解](https://zhuanlan.zhihu.com/p/74736706)
[5] [GAMES-101 图形学入门学习笔记 - 2- 光栅化与着色](https://zhuanlan.zhihu.com/p/140779399)