Chapter 2:Bounding Volume Hierarchies
今天我们来讲层次包围盒, 乍一看比较难, 篇幅也多, 但是咱们一步一步来, 相信大家应该都能听懂
BVH 和 Perlin textures 是这本书中最难的两章, 为什么把 BVH 放在第二章讲呢, 据说, 层次包围盒的渲染效率是比较高的, 考虑到大家的渲染时间开销, 所以先讲这个
光线 - 物体相交是光线跟踪器中的主要时间瓶颈, 运行时间与物体数量成线性关系. 但它是对同一模型的重复搜索, 所以我们可以采用二分搜索的速度进行对数级复杂度搜索. 因为我们在同一模型上发送数百万到数十亿的光线, 我们可以对模型进行分类, 然后每个光线交叉点可以是次线性搜索. 两个最常见的两类是
1)划分空间
2)划分对象
后者通常更容易编码, 并且运行速度与大多数模型一样快. 关键思想是找到一个完全包围 (边界) 所有对象的体积. 例如, 你计算得出了 10 个对象包围盒. 任何未与包围盒相交的射线肯定不会和十个对象相交. 如果射线击中了包围盒, 那么它可能会击中十个物体中的一个. 所以伪代码:
- if (ray hits bounding object)
- return whether ray hits bounded objects
- else
- return false
最关键的一件事情是我们如何将各个物体划分到每个子集中, 单个子集为一个包围盒
引用书上一张图
蓝色和红色边界包含在紫色边界中, 但它们可能重叠, 并且它们不是有序的, 它们只是单纯地被包含在内部.
右边是左图的树结构
对应的检测伪代码是:
- if (hits purple)
- hit0 = hits blue enclosed objects
- hit1 = hits red enclosed objects
- if (hit0 or hit1)
- return true and info of closer hit
- return false
好了, 我们下面就来实现上面的伪代码框架
我们需要将场景中的物体进行划分, 且包围盒需要非常紧凑, 以及考虑光线与包围盒相交的方法, 计算量尽量少. 在大多数模型的实践中, 轴对齐的盒子 (axis - aligned bounding box, 即 AABB) 比较好一些, 我们只需要知道光线是否击中盒体, 而无需操心撞击点的任何信息
有一种常用的 "slab" 的方法, 它是基于 n 个维度的 AABB, 就是取 n 个坐标轴上的区间表示, 称为 "slabs"
一维空间的一个区间, 比如: x∈[3,5] , 它是一条线段
二维空间的一个区间, 比如: x∈[3,5] ,y∈[3,5] , 它是一块矩形区域
我们来确定光线与区间的相交信息
>假设, 上述图中的这种情况, 我们的光线和 x = x0, x = x1 相交于 t0 和 t1
回顾视线方程: p(t) = a+t*b, 若为 x 区间相交, 则方程写为 x(t) = a+ t *b
则 x0 = a.x + t0 * b.x => t0 = (x0 - a.x) / b.x
同理可得 t1 = (x1 - a.x) / b.x
>如果是二维的情况, 那么, 就要加上 y(t) = a+ t *b
我们的工作就是:
计算 Q1(tx0,tx1)
计算 Q2(ty0,ty1)
Q1 和 Q2 是否有交集
大致几种分类情况如下:(下面红色代表 ty, 绿色代表 tx)
综上, 我们发现:
如果两个区间的左端点最大值 小于 右端点的最小值
光线一定和区域有交点
反之
光线和区域相离
>如果是三维的, 那么同理, 步骤如下:
计算 Q1(tx0,tx1)
计算 Q2(ty0,ty1)
计算 Q3(tz0,tz1)
Q1,Q2 和 Q3 是否有交集
代码简单描写为:
_min 和 _max 分别指的是包围盒中三个维度区间的左端点集合和右端点集合
它们是 aabb 的两个数据成员
(这个代码只用来理解用)
后来我们观察发现, 一般情况下满足 t0 <= t1 , 但是有时候 t0> t1
当且仅当视线的方向向量在当前计算维度的分量中为负, 此时 t0> t1
所以, 我们改写成如下形式:
- inline bool aabb::hit(const ray& sight, rtvar tmin, rtvar tmax)const
- {
- for (int i = 0; i <3; ++i)
- {
- rtvar div = 1.0 / sight.direction()[i];
- rtvar t1 = (_min[i] - sight.origin()[i]) / sight.direction()[i];
- rtvar t2 = (_max[i] - sight.origin()[i]) / sight.direction()[i];
- if (div < 0.)stds swap(t1, t2);
- if (stds min(t2, tmax) <= stds max(t1, tmin))
- return false;
- }
- return true;
- }
同时, 我们的相交类也要改一下
我们从现在开始增加各个子类实现
关于 sphere, 球体的三个维度的左端点集合和右端点集合分别为
- aabb sphere::getbox()const
- {
- return aabb(_heart - rtvec(_radius, _radius, _radius), _heart + rtvec(_radius, _radius, _radius));
- }
对于 moving_sphere 我们需要综合开始时刻和结束时刻两个球的盒体的边界
aabb _surrounding_box(aabb box1, aabb box2);
但是, 出于某种考虑, 我觉得把它放在 aabb 盒体类中作为静态成员函数比较好
- /// aabb_box.hpp
- // -----------------------------------------------------
- // [author] lv
- // [begin ] 2019.1
- // [brief ] the aabb-class for the ray-tracing project
- // from the 《ray tracing the next week》
- // -----------------------------------------------------
- namespace rt
- {
- //the statement of aabb class
- class aabb
- {
- public:
- aabb() { }
- aabb(const rtvec& a, const rtvec& b);
- inline bool hit(const ray& sight, rtvar tmin, rtvar tmax)const;
- static aabb _surrounding_box(aabb box1, aabb box2);
- public:
- inline rtvec min()const { return _min; }
- inline rtvec max()const { return _max; }
- private:
- rtvec _min;
- rtvec _max;
- };
- //the implementation of aabb class
- inline aabb::aabb(const rtvec& a, const rtvec& b)
- :_min(a)
- , _max(b)
- {
- }
- inline bool aabb::hit(const ray& sight, rtvar tmin, rtvar tmax)const
- {
- for (int i = 0; i < 3; ++i)
- {
- rtvar div = 1.0 / sight.direction()[i];
- rtvar t1 = (_min[i] - sight.origin()[i]) / sight.direction()[i];
- rtvar t2 = (_max[i] - sight.origin()[i]) / sight.direction()[i];
- if (div < 0.)stds swap(t1, t2);
- if (stds min(t2, tmax) <= stds max(t1, tmin))
- return false;
- }
- return true;
- }
- aabb aabb::_surrounding_box(aabb box1, aabb box2)
- {
- auto fmin = [](const rtvar a, const rtvar b) {return a < b ? a : b; };
- auto fmax = [](const rtvar a, const rtvar b) {return a> b ? a : b; };
- rtvec min{ fmin(box1.min().x(),box2.min().x()),
- fmin(box1.min().y(),box2.min().y()),
- fmin(box1.min().z(),box2.min().z()) };
- rtvec max{ fmax(box1.max().x(),box2.max().x()),
- fmax(box1.max().y(),box2.max().y()),
- fmax(box1.max().z(),box2.max().z()) };
- return aabb(min, max);
- }
- }
- aabb moving_sphere::getbox()const
- {
- rtvec delt{ _radius, _radius, _radius };
- return aabb::_surrounding_box(aabb(_heart1 - delt, _heart1 + delt), aabb(_heart2 - delt, _heart2 + delt));
- }
现在我们开始着手, 划分物体, 并解决 "光线是否击中了当前盒体" 这个开篇的问题
首先, 我们需要创建像开篇那张图中的一颗盒体范围树
树节点定义:
- class bvh_node :public intersect
- {
- public:
- bvh_node() { }
- bvh_node(intersect** world, const int n, const rtvar time1, const rtvar time2);
- virtual bool hit(const ray& sight, rtvar t_min, rtvar t_max, hitInfo& info)const override;
- virtual aabb getbox()const override;
- private:
- intersect* _left;
- intersect* _right;
- aabb _box;
- };
- aabb bvh_node::getbox()const
- {
- return _box;
- }
构造函数中那两个时间实在不知道有什么用(=.=)
之后我们就需要写 hit 函数了, 其实很好写
树结构, 遍历左子树遍历右子树, 返回离 eye 最近的撞击点信息即可
- bool bvh_node::hit(const ray& sight, rtvar t_min, rtvar t_max, hitInfo& info)const
- {
- if (_box.hit(sight, t_min, t_max))
- {
- hitInfo linfo, rinfo;
- bool lhit = _left->hit(sight, t_min, t_max, linfo);
- bool rhit = _right->hit(sight, t_min, t_max, rinfo);
- if (lhit && rhit)
- {
- if (linfo._t <rinfo._t)
- info = linfo;
- else
- info = rinfo;
- return true;
- }
- else if (lhit)
- {
- info = linfo;
- return true;
- }
- else if (rhit)
- {
- info = rinfo;
- return true;
- }
- else
- return false;
- }
- else
- return false;
- }
- }
构造函数设计:
1)随机选择一个轴
2)使用库 qsort 对物体进行排序
3)在每个子树中放一半物体
并且特判了两种情况(物体个数为 1 或者 2)
如果我只有一个元素, 我在每个子树中复制它. 两个物体的话, 一边一个.
明确检查三个元素并且只跟随一个递归可能会有所帮助, 但我认为整个方法将在以后进行优化. 即:
- inline bvh_node::bvh_node(intersect** world, const int n, const rtvar time1, const rtvar time2)
- {
- int axis = static_cast<int>(3 * lvgm::rand01());
- if (axis == 0)
- qsort(world, n, sizeof(intersect*), _x_cmp);
- else if (axis == 1)
- qsort(world, n, sizeof(intersect*), _y_cmp);
- else
- qsort(world, n, sizeof(intersect*), _z_cmp);
- if (n == 1)
- _left = _right = world[0];
- else if (n == 2)
- _left = world[0],
- _right = world[1];
- else
- _left = new bvh_node(world, n / 2, time1, time2),
- _right = new bvh_node(world + n / 2, n - n / 2, time1, time2);
- aabb lbox = _left->getbox();
- aabb rbox = _right->getbox();
- _box = aabb::_surrounding_box(lbox, rbox);
- }
比较函数以_x_cmp 为例:
- inline int _x_cmp(const void * lhs, const void * rhs)
- {
- intersect * lc = *(intersect**)lhs;
- intersect * rc = *(intersect**)rhs;
- aabb lbox = lc->getbox();
- aabb rbox = rc->getbox();
- if (lbox.min().x() - rbox.min().x() <0.)
- return -1;
- else
- return 1;
- }
整个方法在之后可能会优化, 但是目前确实不咋好:
测试代码
- #define LOWPRECISION
- #include <fstream>
- #include "RTmaterial.hpp"
- #include "RThit.hpp"
- #include "camera.hpp"
- using namespace rt;
- int Cnt;
- rtvec lerp(const ray& sight, intersect* world, int depth)
- {
- hitInfo info;
- if (world->hit(sight, (rtvar)0.001, rtInf(), info))
- {
- ray scattered;
- rtvec attenuation;
- if (depth <50 && info.materialp->scatter(sight, info, attenuation, scattered))
- return attenuation * lerp(scattered, world, depth + 1);
- else
- return rtvec(0, 0, 0);
- }
- else
- {
- rtvec unit_dir = sight.direction().ret_unitization();
- rtvar t = 0.5*(unit_dir.y() + 1.);
- return (1. - t)*rtvec(1., 1., 1.) + t*rtvec(0.5, 0.7, 1.0);
- }
- }
- intersect* random_sphere()
- {
- int cnt = 50000;
- intersect **list = new intersect*[cnt + 1];
- list[0] = new sphere(rtvec(0, -1000, 0), 1000, new lambertian(rtvec(0.5, 0.5, 0.5)));
- int size = 1;
- for (int a = -5; a <5; ++a)
- for (int b = -5; b < 5; ++b)
- {
- rtvar choose_mat = lvgm::rand01();
- rtvec center(a + 0.9 * lvgm::rand01(), 0.2, b + 0.9*lvgm::rand01());
- if ((center - rtvec(4, 0.2, 0)).normal()>0.9)
- {
- if (choose_mat <0.55)
- list[size++] = new moving_sphere(center, center + rtvec(0, 0.5*lvgm::rand01(), 0), 0., 1., 0.2,
- new lambertian(rtvec(lvgm::rand01()*lvgm::rand01(), lvgm::rand01()*lvgm::rand01(), lvgm::rand01()*lvgm::rand01())));
- else if (choose_mat < 0.85)
- list[size++] = new sphere(center, 0.2,
- new metal(rtvec(0.5*(1 + lvgm::rand01()), 0.5*(1 + lvgm::rand01()), 0.5*(1 + lvgm::rand01())), 0.5*lvgm::rand01()));
- else
- list[size++] = new sphere(center, 0.2,
- new dielectric(1.5));
- }
- }
- list[size++] = new sphere(rtvec(0, 1, 0), 1.0, new dielectric(1.5));
- list[size++] = new moving_sphere(rtvec(-4.5, 1, 0.65), rtvec(-4,1,0.15), 0., 1., 1.0,
- new lambertian(rtvec(0.4, 0.2, 0.1)));
- list[size++] = new sphere(rtvec(4, 1, 0), 1.0, new metal(rtvec(0.7, 0.6, 0.5), 0.));
- return new intersections(list, size);
- }
- void build_12_1()
- {
- stds ofstream file("graph1-2.ppm");
- size_t W = 200, H = 100, sample = 100;
- if (file.is_open())
- {
- file << "P3\n" << W << "" << H <<"\n255\n" << stds endl;
- intersect* world = random_sphere();
- rtvec lookfrom(13, 2, 3);
- rtvec lookat(0, 0, 0);
- float dist_to_focus = 10.0;
- float aperture = 0.0;
- camera cma(lookfrom, lookat, rtvec(0, 1, 0), 20, rtvar(W) / rtvar(H), aperture, 0.7*dist_to_focus, 0., 1.);
- for (int y = H - 1; y>= 0; --y)
- for (int x = 0; x <W; ++x)
- {
- rtvec color;
- for (int cnt = 0; cnt < sample; ++cnt)
- {
- lvgm::vec2<rtvar> para{
- (lvgm::rand01() + x) / W,
- (lvgm::rand01() + y) / H };
- color += lerp(cma.get_ray(para), world, 0);
- }
- color /= sample;
- color = rtvec(sqrt(color.r()), sqrt(color.g()), sqrt(color.b())); //gamma 校正
- int r = int(255.99 * color.r());
- int g = int(255.99 * color.g());
- int b = int(255.99 * color.b());
- file << r << "" << g <<" " << b << stds endl;
- }
- file.close();
- if (world)delete world;
- stds cout << "complished" << stds endl;
- }
- else
- stds cerr << "open file error" << stds endl;
- }
- int main()
- {
- build_12_1();
- }
- main.cpp
有点伤心且不知所措, 但其实还是学了很多相交的知识的~
感谢您的阅读, 生活愉快~
来源: https://www.cnblogs.com/lv-anchoret/p/10284085.html