Acceleration Structures

我们已经知道我们可以通过AABB来加速我们的光线与物体求交，但是在这个基础上我们可以进一步进行加速，下面是两种方法

Uniform Grid

找到整个场景的AABB包围盒
将整个包围盒均匀划分为网格
找到所有物体在网格中的位置并存储起来

开始进行光线与模型求交，光线与物体求交是非常快的，光线与模型求交是非常慢的。我们的光线会与一些格子相交，根据我们上面在网格中存储的信息，如果我们的光线经过了存储物体信息的格子，那么这个光线就有可能与物体相交，此时我们让光线与这个物体求一次交，如果不存在交点就进行下一个格子。那么我们怎么知道这个光线经过了那么些格子呢？当我们将光线看作一条线段，格子看成我们熟悉的像素，那么这其实就是如何光栅化一条直线，我们可以使用DDA或者Bresenham算法进行判断哪些格子被射线所经过，再查询格子中的物体信息，如果有物体再计算这个物体和光线是否相交。当然我们这是考虑的二维的情况，三维中也有对应的方法。

当然这其中的格子的划分数量也是有一个权衡的，如果格子画的太稀疏很多物体就会在一个格子中，光线就会多次与物体求交，如果格子画的太密集一个物体就会存在在多个格子中，内存占用就会非常大，而且光线与盒子求交的次数也会非常多

在不同的场景中使用均匀网格法的效率也是不同的比如下面这副场景，场景中的物体在场景中分布的比较均匀那么在格子中物体也会被均匀的分配到不同的格子中，那么网格法的效率会更高在下面这个场景来说，场景中间有很多空气，空气中没有物体，物体的分布不均匀，这像是一个运动场中放一个茶壶，但是网格都是茶壶那么大的大小，这样划分的格子就会很多，那么效率就会非常低

Spatial Partitions

很显然均匀网格法存在问题，所以我们需要改进这个方法，我们想既然均匀网格不好，那么我们可不可以在物体密集的地方使用多的网格，在物体少的地方使用少的网格，这样的技术我们称为空间划分。下面是三个不同的空间划分的方法

Oct-Tree 八叉树（二维是四叉树），我们都学过二叉树，那这个方法其实就是类似的，在二维情况下将整个场景使用一个AABB作为根节点，然后在平均分为4个子节点，在子节点中在进行平均划分为4个节点，我们可以定义，如果一个节点中的物体少于n的时候不再进行划分。如图所示，它会存在这种情况：一个物体存在多个节点中
KD-Tree 它是这么做的：每次水平竖直的划分空间交替进行，这样做让他始终是一个二叉树的结构
BSP-Tree 与KD-Tree的区别就是斜着切

KD-Tree

我们的这些空间划分都是预计算操作，也就是在光线追踪前做的我们最常用的就是KD-Tree

建立KD-Tree

我们一步一步看如何建立一个KD-Tree，我们先看看半部分子树的划分（另一部分子树其实也需要同步进行）我们只需要在叶子节点记录物体是否存在的信息，父节点记录节点的关系就行了。当然我们没必要每次在空间的中间进行划分，但我们通常这么做。

如何使用KD-Tree进行光追加速

这个光线首先是和A节点相交那么这个光线有可能和节点1和节点B相交，我们假设我们最终划分的树就是这样，1是树的叶子节点，那么我们判断1的时候因为它是叶子节点，所以我们需要判断光线与叶子节点中的所有物体是否相交因为我们的光线与B是相交的那么光线就有可能与B的子节点2和C相交，我们先判断光线是否与节点2相交，而节点2是叶子节点所以我们又需要让光线与2节点中所有的物体进行求交我们再判断光线是否与C相交，而C不是叶子节点，所以光线有可能与C的子节点相交，那么我们需要判断光线是否与D和节点3相交，依次类推

当然这KD-Tree还存在一个问题，就是你不好判读这个物体是否存在在一个AABB中，因为你需要将物体的id信息存放在叶子节点中，而这就涉及到了一个问题，如何判断一个三角形是否与AABB相交，而这是一个非常耗时，非常困难的操作。所以我们需要寻找新的方法。

Object Partitions

Bounding Volume Hierarchy (BVH)

既然我们从空间的角度划分物体会产生这个问题，那么我们换一个方向，从物体开始划分场景，也就是BVH 这个方法在图形学中得到了非常广泛的应用

我们可以从物体的角度，把所有的物体的AABB作为根节点，让后把物体分成两个部分再求两个部分物体的AABB 同样我们继续，我们可以定义如果一个节点中只有n个物体的时候那么不再进行划分这个方法可以让一个物体只出现在一个叶子节点中，也不用让三角形与AABB求交，当然我们会发现它的AABB会相交，但是好像并没有什么影响

总结一下就是

找到物体的包围盒
将物体分为两个部分（当然划分也是有讲究的）
- 通常我们是水平划分依次，竖直划分一次
- 如果一个场景是在一个长条中，如一个火车，这样的场景我们通常会进行多次竖直划分，可以根据AABB的长宽高来进行选择
- 或者我们可以总是选择最长AABB的长宽高进行划分
- 或者我们可以选择中间的物体进行划分
重新计算出两个部分的包围盒
当划分到当如果一个叶子节点中物体少于n个的时候不再进行划分
将物体的id信息存放在叶子节点中

如何使用BVH进行光追加速

这个过程和KD-Tree类似这里直接给出伪代码

辐射度量学

对于15、17节课我已经有了一篇详细介绍的文章了，这里就不再赘述了

作业5解析

作业要求：

你需要修改的函数是：
Renderer.cpp 中的 Render()：这里你需要为每个像素生成一条对应的光线，然后调用函数 castRay() 来得到颜色，最后将颜色存储在帧缓冲区的相应像素中。
Triangle.hpp 中的 rayTriangleIntersect(): v0, v1, v2 是三角形的三个顶点，orig 是光线的起点，dir 是光线单位化的方向向量。tnear, u, v 是你需要使用我们课上推导的 Moller-Trumbore 算法来更新的参数。

这次的作业非常简单，我们只需要使用框架实现两个函数就行了，而且第二题的公式我们在上一篇文章中已经给出了。

实现后的Render()函数

 1// [comment]
 2// The main render function. This where we iterate over all pixels in the image, generate
 3// primary rays and cast these rays into the scene. The content of the framebuffer is
 4// saved to a file.
 5// [/comment]
 6void Renderer::Render(const Scene& scene)
 7{
 8    std::vector<Vector3f> framebuffer(scene.width * scene.height);
 9
10    float scale = std::tan(deg2rad(scene.fov * 0.5f));
11    float imageAspectRatio = scene.width / (float)scene.height;
12    
13    // Use this variable as the eye position to start your rays.
14    Vector3f eye_pos(0);
15    int m = 0;
16    for (int j = 0; j < scene.height; ++j)
17    {
18        for (int i = 0; i < scene.width; ++i)
19        {
20            // generate primary ray direction
21            float x;
22            float y;
23            // TODO: Find the x and y positions of the current pixel to get the direction
24            // vector that passes through it.
25            // Also, don't forget to multiply both of them with the variable *scale*, and
26            // x (horizontal) variable with the *imageAspectRatio*            
27            x = (2 * (i + 0.5f) / (float)scene.width - 1) * imageAspectRatio * scale;
28            y = (1 - 2 * (j + 0.5f) / (float)scene.height) * scale;
29
30            Vector3f dir = normalize(Vector3f(x, y, -1)); // Don't forget to normalize this direction!
31            framebuffer[m++] = castRay(eye_pos, dir, scene, 0);
32        }
33        UpdateProgress(j / (float)scene.height);
34    }
35
36    // save framebuffer to file
37    FILE* fp = fopen("binary.ppm", "wb");
38    (void)fprintf(fp, "P6\n%d %d\n255\n", scene.width, scene.height);
39    for (auto i = 0; i < scene.height * scene.width; ++i) {
40        static unsigned char color[3];
41        color[0] = (char)(255 * clamp(0, 1, framebuffer[i].x));
42        color[1] = (char)(255 * clamp(0, 1, framebuffer[i].y));
43        color[2] = (char)(255 * clamp(0, 1, framebuffer[i].z));
44        fwrite(color, 1, 3, fp);
45    }
46    fclose(fp);    
47}

我们解释一下为什么x和y要这么算。
首先我们的目的是将屏幕空间的像素坐标一步步转换回相机空间的 3D 坐标，也就是将屏幕空间$[0,0]\times[width,height]$转换成相机空间$[-scale \times imageAspectRatio,-scale]\times[scale \times imageAspectRatio,scale]$，这里相机空间的y为什么不用乘以imageAspectRatio呢？因为这是规定好的，我们使用y去定义了视野的上下边界。

对x
- 首先我们在屏幕空间中$\frac{i + 0.5}{width}$ 将坐标映射到 $[0, 1]$
- 乘以 2 再减 1，即 $2 \cdot [0, 1] - 1$，将范围变换到 $[-1, 1]$
- 乘以 scale，即 $scale \cdot [0, 1] - 1$，将范围变换到 $[-scale, scale]$
- 最后再乘以 imageAspectRatio，即 $scale \cdot [0, 1] - 1$，将范围变换到 $[-scale \times imageAspectRatio, scale \times imageAspectRatio]$
对y
- 在屏幕空间中，坐标原点 $(0,0)$ 通常在左上角，$j$ 增大是向下的，为了让屏幕上方的像素对应 3D 空间中正的 $y$ 值，我们需要进行反转
- 我们首先和x一样，将坐标映射到 $[0, 1]$ 使用$\frac{j + 0.5}{height}$
- 然后我们对它乘以2 变成了 $[0, 2]$ 在用1减去y 得到$[1, -1]$，我们就对它进行了反转了
- 然后乘以scale，即 $scale \cdot [1, -1]$，将范围变换到 $[-scale, scale]$

总结公式

$$x = \underbrace{\left( \frac{2 \cdot (i + 0.5)}{width} - 1 \right)}_{\text{Range: } [-1, 1]} \cdot \text{aspectRatio} \cdot scale$$

$$y = \underbrace{\left( 1 - \frac{2 \cdot (j + 0.5)}{height} \right)}_{\text{Range: } [1, -1]} \cdot scale$$

最后对向量归一化就行了

实现后的rayTriangleIntersect()函数根据我们之前的公式

 1
 2bool rayTriangleIntersect(const Vector3f &v0, const Vector3f &v1,
 3                          const Vector3f &v2, const Vector3f &orig,
 4                          const Vector3f &dir, float &tnear, float &u,
 5                          float &v) {
 6  // TODO: Implement this function that tests whether the triangle
 7  // that's specified bt v0, v1 and v2 intersects with the ray (whose
 8  // origin is *orig* and direction is *dir*)
 9  // Also don't forget to update tnear, u and v.
10  constexpr float kEpsilon = 1e-8f;
11  const Vector3f e1 = v1 - v0;
12  const Vector3f e2 = v2 - v0;
13  const Vector3f s1 = crossProduct(dir, e2);
14  const float det = dotProduct(e1, s1);
15
16  // 判断射线和平面是否平行
17  if (std::fabs(det) < kEpsilon) {
18    return false;
19  }
20
21  const float invDet = 1.0f / det;
22  const Vector3f s = orig - v0;
23  u = dotProduct(s, s1) * invDet;
24
25  // 约束重心坐标u在[0, 1]范围内
26  if (u < 0.0f || u > 1.0f) {
27    return false;
28  }
29
30  const Vector3f s2 = crossProduct(s, e1);
31  v = dotProduct(dir, s2) * invDet;
32
33  // 约束重心坐标v在[0, 1]范围内，并且u + v <= 1，保证交点在三角形内部
34  if (v < 0.0f || u + v > 1.0f) {
35    return false;
36  }
37
38  tnear = dotProduct(e2, s2) * invDet;
39  return tnear > kEpsilon;
40}

照着抄就行
最后的到的效果是

Acceleration Structures#

Uniform Grid#

Spatial Partitions#

KD-Tree#

建立KD-Tree#

如何使用KD-Tree进行光追加速#

Object Partitions#

Bounding Volume Hierarchy (BVH)#

如何使用BVH进行光追加速#

辐射度量学#

对于15、17节课我已经有了一篇详细介绍的文章了，这里就不再赘述了#

作业5解析#