前言
- 本系列总结自Ray Tracing: The Next Week,前篇为rayTracingInOneWeekend – 爱莉希雅 – 博客园 (cnblogs.com)。需要一点点的图形学基础,不懂的同学可以移至GAMES101-现代计算机图形学入门-闫令琪_哔哩哔哩_bilibili
- 这个版本的光线追踪非常适合初学者,难度既不大但也感受到光追的魅力
总结
动态模糊
动态模糊:在生活中,摄像机快开打开的时间间隔中,摄像机和物体都在移动,拍出来的图像肯定是糊的。大多数特效都可以通过对每个像素进行多重采样实现,动态模糊也属于这个范畴。
方法:我们可以
随机
的方法在不同时间发射多条射线来模拟快门的打开
。也就是说,在快门打开时,随时间变化随机生成光线,并同时射出与模型相交。一般而言,我们让摄像机和物体同时运动,并让每条射线都拥有自己的时间点
只要物体在那个时间处于它自己正确的位置,我们便能得出光线在那个时间的精确平均值
实现。让漫反射材质的球动起来:
让每条射线拥有自己的时间点
class ray{ public: ray() {} ray(const vec3& origin, const vec3& direction, double time = 0.0) : orig(origin), dir(direction), tm( time ) {} vec3 get_origin() const { return orig; } vec3 get_direction() const { return dir; } vec3 at(double t) const { return orig + dir * t; } private: vec3 orig; vec3 dir; double tm;};
让camera在[time1, time2]中随机生成光线
class camera{ public: camera(vec3 lookfrom, vec3 lookat, vec3 up_vector, double boc, double aspect, double aperture, double focus_dist, double t0 = 0.0, double t1 = 0.0) { origin = lookfrom; lens_radius = aperture / 2; time0 = t0; time1 = t1; auto theta = degrees_to_radians(boc); auto half_height = tan(theta / 2); auto half_width = aspect * half_height; w = unit_vector(lookfrom - lookat); u = unit_vector(cross(up_vector, w)); v = cross(w, u); lower_left_corner = origin - u * half_width * focus_dist - v * half_width * focus_dist - w * focus_dist; horizontal = u * 2 * half_width * focus_dist; vertical = v * 2 * half_height * focus_dist; } ray get_ray(double s, double t) { vec3 rd = random_in_unit_disk() * lens_radius; vec3 offset = u * rd.x() + v * rd.y(); return ray(origin + offset, lower_left_corner + horizontal * s + vertical * t - origin - offset, random_double(time0, time1)); } private: vec3 lower_left_corner; vec3 horizontal; vec3 vertical; vec3 origin; vec3 u, v, w; double lens_radius; double time0, time1;};
让sphere class的球心在[time0, time1]内从center0运动至center1。超出这个时间范围,球心依然在动,也就是说,线性插值time既可小于0.0,亦可大于1.0。这两个时间变量和快门的开关时刻无需一一对应
class moving_sphere : public hittable{ public: moving_sphere() = default; moving_sphere(vec3 cen0, vec3 cen1, double t0, double t1, double r, shared_ptr m_ptr) : center0(cen0), center1(cen1), time0(t0), time1(t1), radius(r), mat_ptr(m_ptr) {} virtual bool hit(const ray& r, const double t_min, const double t_max, hit_record& rec) const; vec3 center(double time) const; private: vec3 center0, center1; double time0, time1; double radius; shared_ptrmat_ptr;};inline vec3 moving_sphere::center(double time) const{ //线性插值 return center0 + ( (center1 - center0) * (time - time0) / (time1 - time0) );}//与sphere的变化只在centerinline bool moving_sphere::hit(const ray& r, const double t_min, const double t_max, hit_record& rec) const{ vec3 oc = r.get_origin() - center(r.get_time()); auto a = r.get_direction().length_squared(); auto half_b = dot(r.get_direction(), oc); auto c = oc.length_squared() - radius * radius; auto delta = half_b * half_b - a * c; if (delta > 0) { auto root = sqrt(delta); auto temp = (-half_b - root) / a; if (temp > t_min && temp t_min && temp < t_max) { rec.t = temp; rec.p = r.at(rec.t); vec3 outward_normal = (rec.p - center(r.get_time())) / radius; rec.set_face_normal(r, outward_normal); rec.mat_ptr = mat_ptr; return true; } } return false;}
确保moving_sphere的材质在运算散射时,散射光线与入射光线存在的时间点相同
//漫反射材质class lambertian : public material{ public: lambertian( const vec3& a) : albedo(a) {} virtual bool scatter(const ray& r_in, const hit_record& rec, vec3& attenuation, ray& scattered) const { vec3 scatter_direction = rec.normal + random_unit_vector(); scattered = ray(rec.p, scatter_direction, r_in.get_time()); attenuation = albedo; return true; } private: vec3 albedo;};
快门在time0打开,在time1关闭,每个球心从位置L运动到L + (0,r/2, 0),r为[0,1]的随机数。跳动肯定时变换y轴
hittable_list random_scene(){ hittable_list world; world.add(make_shared(vec3(0, -1000, 0), 1000, make_shared(vec3(0.5, 0.5, 0.5)))); int i = 1; for (int a = -10; a < 10; ++a) { for (int b = -10; b 0.9) { if (choose_mat < 0.8) { //diffuse auto albedo = vec3::random() * vec3::random(); world.add(make_shared(center, center + vec3(0, random_double(0, 0.5), 0), 0.0, 1.0, 0.2, make_shared(albedo))); } else if (choose_mat < 0.95) { //metal auto albedo = vec3::random(0.5, 1); auto fuzz = random_double(0, 0.5); world.add(make_shared(center, 0.2, make_shared(albedo, fuzz))); } else { //glass world.add(make_shared(center, 0.2, make_shared(1.5))); } } } } world.add(make_shared(vec3(0, 1, 0), 1.0, make_shared(1.5))); world.add( make_shared(vec3(-4, 1, 0), 1.0, make_shared(vec3(0.4, 0.2, 0.1)))); world.add( make_shared(vec3(4, 1, 0), 1.0, make_shared(vec3(0.7, 0.6, 0.5), 0.0))); return world;}void output_image(){ //... auto world = random_scene(); vec3 lookfrom(13, 2, 3); vec3 lookat(0, 0, 0); vec3 up_vector(0, 1, 0); auto dist_to_focus = 10.0; auto aperture = 0.0; const auto aspect_ratio = double(image_width) / image_height; //...}
输出:
总结:用
随机
的方法在不同时间发射多条射线来模拟快门的打开
。让摄像机和物体同时运动,并让每条射线都拥有自己的时间点
层次包围盒
在这之前,光线的求交运算会耗费大量时间,且计算所消耗的时间与场景中的物体数量线性相关。使用遍历反复查找同一个模型会有许多多余的计算,因此本节将使用一种
树型结构
来进行加速,即层次包围盒
(Bounding Volume Hierarchy)层次包围盒又属于
空间数据结构
,空间数据结构有以下优点
:- 易于将二维和三维的概念扩展到高维
- 可用于加速查询。如裁减算法、相交测试、光线追踪、碰撞检测等。访问速度很快。时间复杂度可以从 O(n)提升到 O(log n)
层次包围盒的思想:用
体积略大而几何特征简单
的包围盒来近似地描述复杂的几何对象,从而只需对包围盒重叠的对象进行进一步的相交测试
。此外,通过构造树状层次结构,可以越来越逼近对象的几何模型,直到几乎完全获得对象的几何特征。不再以空间作为划分依据,而是从对象的角度考虑
层次包围盒的代码大致如下:
if (ray hits bounding object) return whether ray hits bounded objectselse return false
层次包围盒的缺点:对于复杂的模型,效果并不好
举个例子,以下场景以层次树结构进行组织,包含一个根节点(root)、一些内部节点(internal nodes),以及一些叶子节点 (leaves)。顶部的节点是根,其无父节点。叶子节点(leaf node)包含需渲染的实际几何体,且其没有子节点
AABB包围盒求交
运用层次包围盒的前提是 如何判断光线是否与物体相交,这里就用到了
AABB包围盒
:由三对平面的交集构成,且任意一对平面都与x轴,y轴、z轴垂直如图:
思想
: 光线若不会与一个物体相交,那就没有必要去遍历该物体的所有三角形面。因此利用一个包围盒包住该物体,在与该物体计算求交之前先判断光线是否与包围盒相交,若包围盒与光线没有交点,那么显然不会与物体有交点举个例子:我们以2维的AABB包围盒为例
如下图,左边的图计算
x轴
上 最晚进来的交点t\(min\) 和 最先出去的t\(max\);中间的图计算y轴
上 最先进来的交点t\(min\)和最先出去的t\(max\);右边的图则综合前两张图的计算结果,计算出最终的结果。最终结果是如何计算的呢?主要是以下几点:- 只有当光线进入了所有平面,才是真正地进入了包围盒。\(\large t_{enter} = max(tmin)\)
- 只有当光线离开了任一平面,就算是真正地离开了包围盒。\(\large t_{exit} = min(t_{max})\)
- 当且仅当,\(\large t_{enter} 0\)成立时,才有交点。
- 若t小于0,说明盒子在光线背后;若\(\large t_{exit} >= 0, t_{enter} < 0\),说明光线起点在盒子内部
为什么使用AABB?因为求t时,只需使用某一轴的信息,不必使用整个坐标,与点乘相比更容易
伪代码实现:
compute(tx0, tx1)compute(ty0, ty1)compute(tz0, tz1)return overlap?( (tx0, tx1), (ty0, ty1), (tz0, tz1))
实现:
若光线向x轴负方向射出并与包围盒相交,很显然此时\(t_{enter} > t_{exit}\),因此我们需要交换一下
而对于Bx = 0来说,此时答案是无意义的,为NaN,因此对于任何NaN结果来说,我们都返回false;对于(x0 – Ax)和(x1 – Ax)两分子为0来说,此时只有一个解,相当于是光线擦边,认定光线射中或没射中都可以,我们认定它没有射中,即return false
class aabb{ private: vec3 m_min; vec3 m_max; public: aabb() = default; aabb(const vec3& a, const vec3& b) { m_min = a; m_max = b; } vec3 get_min() const { return m_min; } vec3 get_max() const { return m_max; } bool hit(const ray& r, double tmin, double tmax) const { for (int i = 0; i < 3; ++i) { //若光线向x轴负方向射出 //库函数考虑了异常、NaN,速度要慢些 //以下考虑了分母为0的情况 auto t0 = ffmin((m_min[i] - r.get_origin()[i]) / r.get_direction()[i], (m_max[i] - r.get_origin()[i]) / r.get_direction()[i]); auto t1 = ffmax((m_min[i] - r.get_origin()[i]) / r.get_direction()[i], (m_max[i] - r.get_origin()[i]) / r.get_direction()[i]); //只有当光线进入了所有平面,才是真正地进入了包围盒 tmin = ffmax(t0, tmin); //只有当光线离开了任一平面,就算是真正地离开了包围盒 tmax = ffmin(t1, tmax); if (tmax <= tmin) { return false; } } return true; }};
每个图元构建aabb的方式。在hittable class中加入bounding_box()来计算包裹hittable class的包围盒,再构建一个层次包围盒的树形结构。在层次树中,所有
图元
都处于叶子节点
,又因为物体会动,所以还需接受time1和time2class hittable{ public: virtual bool hit(const ray& r, const double t_min, const double t_max, hit_record& rec) const = 0; virtual bool bounding_box(double t0, double t1, aabb& output_box) const = 0;};//计算两个包围盒的包围盒aabb surrounding_box(aabb box0, aabb box1){ vec3 small(ffmin(box0.get_min().x(), box1.get_min().x()), ffmin(box0.get_min().y(), box1.get_min().y()), ffmin(box0.get_min().z(), box1.get_min().z())); vec3 big(ffmax(box0.get_max().x(), box1.get_max().x()), ffmax(box0.get_max().y(), box1.get_max().y()), ffmax(box0.get_max().z(), box1.get_max().z())); return aabb(small, big);}//sphere构建包围盒很简单bool sphere::bounding_box(double t0, double t1, aabb& output_box) const{ output_box = aabb(center - vec3(radius, radius, radius), center + vec3(radius, radius, radius)); return true;}//moving_sphere构建包围盒,先求球体在t0时的包围盒,再求t1时的包围盒,再计算这两个包围盒的包围盒bool moving_sphere::bounding_box(double t0, double t1, aabb& output_box) const{ aabb box0(center(t0) - vec3(radius, radius, radius), center(t0) + vec3(radius, radius, radius)); aabb box1(center(t1) - vec3(radius, radius, radius), center(t1) + vec3(radius, radius, radius)); output_box = surrounding_box(box0, box1); return true;}//hittable_list构建包围盒,在运行时计算(构造函数中进行运算亦可),因为包围盒的计算一般只有在BVH构造时才调用inline bool hittable_list::bounding_box(double t0, double t1, aabb& output_box) const{ if (objects.empty() == true) { return false; } aabb temp_box; bool first_box = true; for (const auto& object : objects) { //图元无法构建包围盒 if (object->bounding_box(t0, t1, temp_box) == false) { return false; } //是否是第一个包围盒,不是则进行构造两个包围盒的包围盒 output_box = first_box ? temp_box : surrounding_box(output_box, temp_box); first_box = false; } return true;}
构建层次包围盒
构建层次包围盒的
分割
原则。每次分割时沿轴把物体列表分成两半:- 随机选取任一一轴进行分割
- 使用sort()对图元进行排序,也就是对图元任一轴的坐标大小进行排序
- 根据排序对半分,每个子树分一半的图元
BVH也应为hittable一员。两种设计方案,我选择前一种:
为树和树的节点设计两个不同的class
一个class加上指针
//constantAndTool.hinline int random_int(int min, int max){ return static_cast(random_double(min, max + 1));}//hittable.h//设计bvh树节点class bvh_node : public hittable{ private: shared_ptrleft; shared_ptrright; aabb box; public: bvh_node() = default; bvh_node(hittable_list& list, double time0, double time1) : bvh_node(list.get_objects(), 0, list.get_objects().size(), time0, time1) {} bvh_node(vector<shared_ptr>& objects, size_t start, size_t end, double time0, double time1); virtual bool hit(const ray& r, const double t_min, const double t_max, hit_record& rec) const; virtual bool bounding_box(double t0, double t1, aabb& output_box) const;};//若光线打中包围盒,则继续判断的子节点inline bool bvh_node::hit(const ray& r, const double t_min, const double t_max, hit_record& rec) const{ //终止条件 if (box.hit(r, t_min, t_max) == false) { return false; } //左右中 bool hit_left = left->hit(r, t_min, t_max, rec); bool hit_right = right->hit(r, t_min, hit_left == true ? rec.t : t_max, rec); return hit_left || hit_right;}//std::sort的compare函数//先判断是哪个轴,再为compare函数赋值inline bool box_compare(const shared_ptr a, const shared_ptrb, int axis){ aabb box_a; aabb box_b; //两图元中任一图元无法构建包围盒 if (a->bounding_box(0, 0, box_a) == false || b->bounding_box(0, 0, box_b) == false) { cerr << "No bounding box in bvh_node constructor." << endl; } return box_a.get_min().e[axis] < box_b.get_min().e[axis];}inline bool box_x_compare(const shared_ptra, const shared_ptr b){ return box_compare(a, b, 0);}inline bool box_y_compare(const shared_ptra, const shared_ptr b){ return box_compare(a, b, 1);}inline bool box_z_compare(const shared_ptra, const shared_ptr b){ return box_compare(a, b, 2);}//构建bvh//当数组传入的元素个数只剩下两个时,在左右子树各放一个;当只有一个元素时,为了省去检查空指针这步,在左右子树都放一个inline bvh_node::bvh_node(const vector<shared_ptr>& objects, size_t start, size_t end, double time0, double time1){ auto temp_objects = objects; int axis = random_int(0, 2); auto comparator = (axis == 0) ? box_x_compare : (axis == 1) ? box_y_compare : box_z_compare; //当前还需分割多少图元 size_t object_span = end - start; if (object_span == 1) { left = right = temp_objects[start]; } else if (object_span == 2) { if (comparator(temp_objects[start], temp_objects[start + 1]) == true) { left = temp_objects[start]; right = temp_objects[start + 1]; } else { left = temp_objects[start + 1]; right = temp_objects[start]; } } else { std::sort(temp_objects.begin() + start, temp_objects.begin() + end, comparator); auto mid = start + object_span / 2; left = make_shared(temp_objects, start, mid, time0, time1); right = make_shared(temp_objects, mid, end, time0, time1); } aabb box_left, box_right; //检验图元是否有aabb if (left->bounding_box(time0, time1, box_left) == false || right->bounding_box(time0, time1, box_right) == false) { cerr << "No bounding box in bvh_node constructor." << endl; } box = surrounding_box(box_left, box_right);}
在random_scene()中进行运用
hittable_list random_scene(){ //... return static_cast(make_shared(world, 0.0, 1.0));}
输出:结果图和上一小节一样但是速度更快
总结:
纹理
纹理贴图:通过将投影方程(projector function)运用于空间中的点 ,从而得到一组称为参数空间值(parameterspacevalues)的关于纹理的数值
纹理贴图常常意味着一个将颜色赋予图元表面的过程,这个过程可以是纹理生成代码,亦可以是一张图片,亦可以是前两者的结合
棋盘格上的球体
texture class。两种设计思路:
- 把静态rgb颜色和贴图写成两
种
class,以此来区分 - 把静态rgb颜色和贴图写成同一
类别
class,以此将任何颜色弄成纹理
class texture{ public: virtual vec3 value(double u, double v, const vec3& p) const = 0;};class constant_texture : public texture{ public: constant_texture() = default; constant_texture( vec3 c ) : color(c) {} virtual vec3 value(double u, double v, const vec3& p) const override { return color; } private: vec3 color;};
- 把静态rgb颜色和贴图写成两
更新hit_record,让其存储相交点的uv信息
struct hit_record{ //... double t; double u; double v; bool front_face;//法线朝向//...};
更新漫反射材质,对其进行纹理贴图
将vec3改为纹理指针
class lambertian : public material{ public: lambertian( shared_ptr a ) : albedo(a) {} virtual bool scatter(const ray& r_in, const hit_record& rec, vec3& attenuation, ray& scattered) const { vec3 scatter_direction = rec.normal + random_unit_vector(); scattered = ray(rec.p, scatter_direction, r_in.get_time()); attenuation = albedo->value( rec.u, rec.v, rec.p); return true; } private: shared_ptr albedo;};
漫反射材质呈现棋盘格纹理
运用sine和cosine函数周期性的变化,在三个维度都乘上周期函数
class checker_texture : public texture{ public: checker_texture() = default; checker_texture( shared_ptr t0, shared_ptr t1 ) : even(t0), odd(t1){} //use the sine and cosine functions to change periodically to create a checkerboard texture virtual vec3 value(double u, double v, const vec3& p) const override { auto sine = sin(10 * p.x()) * sin(10 * p.y()) * sin(10 * p.z()); if (sine value(u, v, p); } else { return even->value(u, v, p); } } private: //The odd and even pointer points to a static texture shared_ptr odd; shared_ptr even;};
更新random_scene()
hittable_list random_scene(){ //... auto checkerboard = make_shared(make_shared(vec3(0.2, 0.3, 0.1)),make_shared(vec3(0.9, 0.9, 0.9))); world.add(make_shared(vec3(0, -1000, 0), 1000, make_shared(checkerboard))); //... if (choose_mat < 0.8) { //diffuse auto albedo = vec3::random() * vec3::random(); world.add(make_shared(center, center + vec3(0, random_double(0, 0.5), 0), 0.0, 1.0, 0.2, make_shared(make_shared(albedo)))); } //... world.add(make_shared(vec3(-4, 1, 0), 1.0, make_shared(make_shared(vec3(0.4, 0.2, 0.1)))));}
输出:
棋盘格球体
更换场景
hittable_list two_spheres(){ hittable_list objects; auto checkerboard = make_shared ( make_shared(vec3(0.2, 0.3, 0.1)), make_shared(vec3(0.9, 0.9, 0.9)) ); objects.add(make_shared(vec3(0, -10, 0), 10, make_shared(checkerboard))); objects.add(make_shared(vec3(0, 10, 0), 10, make_shared(checkerboard))); return static_cast(make_shared(objects, 0.0, 1.0));}void output_image(){ //... auto world = two_spheres(); vec3 lookfrom(13, 2, 3); vec3 lookat(0, 0, 0); vec3 up_vector(0, 1, 0); auto dist_to_focus = 10.0; auto aperture = 0.0; const auto aspect_ratio = double(image_width) / image_height; camera cam(lookfrom, lookat, up_vector, 20, aspect_ratio, aperture, dist_to_focus, 0.0, 1.0);}
输出:
总结
:让材质接受纹理(vec3),根据记录的纹理参数(hit_record)进行散射运算
柏林噪声概念
什么是柏林噪声?
大多数
程序贴图纹理
(计算机算法生成的)都是基于噪声函数的,而柏林噪声就属于噪声函数,这些程序贴图纹理旨在创建用于纹理映射的自然元素(大自然中的一些纹理)的真实表面或三维物体创建的纹理图像为什么使用柏林噪声?
普通噪声毫无规律可言,而想要模拟一些
逼真的随机效果
需要在随机中又有规律以下是普通噪声:
柏林噪声可以生成伪随机数,但其并不是完全随机,对其加以一些物体的规律的限制,可以用于模拟大理石的纹理、云、火焰等(这些自然现象都是随机的,但随机中又有规律)
噪声函数的目的:在三维空间中,提供一种高效地实现、
可重复
、伪随机的信号,此信号大部分能量集中在一个空间频率附近,而视觉上各向同性的(统计上不随旋转变化)以下是柏林噪声:创建某种类似于完全的随机信号,通过低通滤波后,滤除所有的空间高频率而变得
模糊
柏林噪声的特点
- 所有视角细节都是相同的大小
- 可复现性。对于每个输入位置,提供一个可重复的伪随机值
- 实现简单快捷
- 有一个已知的范围(通常是[-1,1])
- 空间频率带宽有限
- 对于重复表现的不明显
- 空间频率在平移下是不变的
如何实现柏林噪声?
三个步骤:
初始化
在3维空间的每个整数(x,y,z)位置产生一个可重复的伪随机值,这运用哈希函数实现,且哈希函数基于一个包含随机数在[0,255]的整数排列表
排列表(Permutation Table):
乱序
存放一系列 索引梯度表(Gradient Table):存放一系列随机梯度值
建立采样空间,对于噪声图上的像素,找到它在单位矩形/单位立方体上对应的一点,并求出它的参考点的坐标
- 对于一维柏林噪声,采样空间为一个一维的坐标轴,轴上整数坐标位置均有一个点,参考点为两侧最近的整数点
- 对于二维柏林噪声,采样空间为一个二维坐标系,坐标系中横纵坐标为整数的地方均有一点,参考点为组成包围该点的单位正方体的四个点
- 三维同理,参考点为组成包围该点的单位立方体的八个点
对于不同类型的噪声,采样点在不同空间中, 根据伪随机整数的梯度来计算 采样空间中小数位置到参考点的距离向量和梯度向量的点积,最后进行插值
一维噪声需一次插值,二维需三次插值,三维需七次插值,呈指数增长,时间复杂度是\(O(2^N)\)
一维
柏林噪声初始化
决定一个点的梯度,我们需要结合
哈希函数
,以点坐标作为参数,以所得值作为索引,在排列表中取得对应值。也就是说,将点坐标转换为排列表的索引,最后取得索引对应的梯度表的值//以下是二维的哈希函数int perm = {...};int grad = {...};void hash( int[] gradient, int x, int y){ //找到对应值在排列表中的索引 int permIdx[] = new int[2]; permIdx[0] = FindIndex(x); permIdx[1] = FindIndex(y); //通过排列表的索引值查找梯度表的索引 int gradIdx[] = new int[2]; gradIdx[0] = perm[permIdx[0]]; gradIdx[0] = perm[permIdx[1]]; //通过梯度表的索引找到对应的梯度值 gradient[0] = grad[gradIdx[0]]; gradient[1] = grad[gradIdx[1]];}
在一维中,各个整点的梯度可以用斜率来表现,因此我们可以不使用梯度表,直接将排列表对应索引取得的值作为梯度
在经典噪声算法中,使用的排列表是由数字0到255散列组成的
int permutation[] = { 151, 160, 137, 91, 90, 15, 131, 13, 201, 95, 96, 53, 194, 233, 7, 225, 140, 36, 103, 30, 69, 142, 8, 99, 37, 240, 21, 10, 23, 190, 6, 148, 247, 120, 234, 75, 0, 26, 197, 62, 94, 252, 219, 203, 117, 35, 11, 32, 57, 177, 33, 88, 237, 149, 56, 87, 174, 20, 125, 136, 171, 168, 68, 175, 74, 165, 71, 134, 139, 48, 27, 166, 77, 146, 158, 231, 83, 111, 229, 122, 60, 211, 133, 230, 220, 105, 92, 41, 55, 46, 245, 40, 244, 102, 143, 54, 65, 25, 63, 161, 1, 216, 80, 73, 209, 76, 132, 187, 208, 89, 18, 169, 200, 196, 135, 130, 116, 188, 159, 86, 164, 100, 109, 198, 173, 186, 3, 64, 52, 217, 226, 250, 124, 123, 5, 202, 38, 147, 118, 126, 255, 82, 85, 212, 207, 206, 59, 227, 47, 16, 58, 17, 182, 189, 28, 42, 223, 183, 170, 213, 119, 248, 152, 2, 44, 154, 163, 70, 221, 153, 101, 155, 167, 43,172, 9, 129, 22, 39, 253, 19, 98, 108, 110, 79, 113, 224, 232, 178, 185, 112, 104, 218, 246, 97, 228, 251, 34, 242, 193, 238, 210, 144, 12, 191, 179, 162, 241, 81, 51, 145, 235, 249, 14, 239, 107, 49, 192, 214, 31, 181, 199, 106, 157, 184, 84, 204, 176, 115, 121, 50, 45, 127, 4, 150, 254, 138, 236, 205, 93, 222, 114, 67, 29, 24, 72, 243, 141, 128, 195, 78, 66, 215, 61, 156, 180 };
采样、插值、计算梯度向量和其偏移量之间的点积
采样:一维柏林噪声函数的参数 输入点 是坐标轴上的一点,其参考点就是它
左右两侧最近
的整数点。因此,我们需要找到这两个点,并求出这两个点对应的梯度
及相对输入点的方向向量
插值:对
梯度值
进行插值,公式为\(\large x = x1 + t(x2 – x1)\),不如线性插值效果不太好,产生的噪声曲线很尖锐。我们使用\(\large t = 3t^2 – 2t^3\)来替换线性插值公式里的\(t\)float perlin( float x ){ //假设x1为x左边的索引值且x - x1的方向向量为正,x2为右边的索引值 int x1 = floor(x); int x2 = x1 + 1; //方向向量 float vec1 = x - x1; float vec2 = x - x2; //求梯度值 float grad1 = perm[x1 % 255] * 2.0 - 255.0; float grad2 = perm[x2 % 255] * 2.0 - 255.0; //插值 float t = 3 * pow( vec1, 2 ) - 2 * pow(vec1, 3); float product1 = grad1 * vec1; float product2 = grad2 * vec2; return product1 + t * (product2 - product1);}
三维柏林噪声
- 将空间(排列表)中的每一点(i,j,k)(坐标为整数)赋值为0,并赋予可以通过
散列
得到的伪随机斜率 - 定义坐标(x,y,z)为一个整数加上小数部分,\((x,y,z) = (i + u, j + v, k + w)\),找到围绕这个点的
单位
立方体的8个参考点,\((i,j,k),(i+1,j,k),·····(i+1, j+1, k+1)\) - 采样、七次线性插值
如下图所示:
- 将空间(排列表)中的每一点(i,j,k)(坐标为整数)赋值为0,并赋予可以通过
柏林噪声的缺点:
- 插值函数每个维度都是\(\large t = 3t^2 – 2t^3\),其二次导数含有非0值。如此在使用噪声函数的导数会导致失真的视觉效果
- 从(i,j,k)得到斜率,需要在一组256个伪随机斜率矢量中查找,这组矢量散布在3D立体球面上。这种分布上的不规则,在噪声函数的结果中产生了不希望的高频
实现大理石材质
将柏林噪声运用在光线追踪中
用一个随机生成的三维数组铺满整个空间,会得到重复的区块
实现:
用哈希表完成
贴瓷砖
class noise_texture : public texture{ private: perlin noise; public: noise_texture() = default; virtual vec3 value(double u, double v, const vec3& p) const override { return vec3(1, 1, 1) * noise.noise(p); }};//perlin noiseclass perlin{ private: static const int point_count = 256; double* ranfloat;//Gradients table //Permutation Table int* perm_x; int* perm_y; int* perm_z; //make the value of permutation table out of order static void permute(int* p, int n) { for (int i = n - 1; i > 0; --i) { int target = random_int(0, i); int tmp = p[i]; p[i] = p[target];//得到随机值 p[target] = tmp;//防止排列表重复 } } //生成排列表对应的伪随机斜率 static int* perlin_generate_perm() { auto p = new int[point_count]; for (int i = 0; i < perlin::point_count; ++i) { p[i] = i; } permute(p, point_count); return p; } public: perlin() { ranfloat = new double[point_count]; for (int i = 0; i < point_count; ++i) { ranfloat[i] = random_double(); } perm_x = perlin_generate_perm(); perm_y = perlin_generate_perm(); perm_z = perlin_generate_perm(); } ~perlin() { delete[] ranfloat; delete[] perm_x; delete[] perm_y; delete[] perm_z; } //算得位置p的伪随机值 double noise(const vec3& p) const { //小数部分 auto u = p.x() - floor(p.x()); auto v = p.y() - floor(p.y()); auto w = p.z() - floor(p.z()); //整数部分 //限制索引值在[0,255] auto i = static_cast(4 * p.x()) & 255; auto j = static_cast(4 * p.y()) & 255; auto k = static_cast(4 * p.z()) & 255; return ranfloat[perm_x[i] ^ perm_y[j] ^ perm_z[k]]; }};hittable_list two_perlin_spheres(){ hittable_list objects; auto pertext = make_shared(); objects.add(make_shared(vec3(0, -1000, 0), 1000, make_shared(pertext))); objects.add(make_shared(vec3(0, 2, 0), 2, make_shared(pertext))); return static_cast(make_shared(objects, 0.0, 1.0));}
- 输出:
我们可以让其更加平滑去掉格子,使用
线性插值
三线性插值:
//三线性插值inline double trilinear_interp(double c[2][2][2], double u, double v, double w){ auto lerp_x_y0_z0 = c[0][0][0] + (u - 0) / 1 * (c[1][0][0] - c[0][0][0]); auto lerp_x_y1_z0 = c[0][1][0] + u * (c[1][1][0] - c[0][1][0]); auto lerp_x_y0_z1 = c[0][0][1] + u * (c[1][0][1] - c[0][0][1]); auto lerp_x_y1_z1 = c[0][1][1] + u * (c[1][1][1] - c[0][1][1]); auto lerp_y_z0 = lerp_x_y0_z0 + v * (lerp_x_y1_z0 - lerp_x_y0_z0); auto lerp_y_z1 = lerp_x_y0_z1 + v * (lerp_x_y1_z1 - lerp_x_y0_z1); auto lerp_z = lerp_y_z0 + w * (lerp_y_z1 - lerp_y_z0); return lerp_z;}class perlin{ //... //Calculate the random value of position P double noise(const vec3& p) const { //decimal part auto u = p.x() - floor(p.x()); auto v = p.y() - floor(p.y()); auto w = p.z() - floor(p.z()); //integral part int i = floor(p.x()); int j = floor(p.y()); int k = floor(p.z()); double c[2][2][2]; //xor //Record the Hasche transform to get the gradient of the reference point //The limiting gradient value is [0,255] for (int di = 0; di < 2; ++di) { for (int dj = 0; dj < 2; ++dj) { for (int dk = 0; dk < 2; ++dk) { //hash transform //Finds the index of the gradient table bythe index values of the permutation table //The gradient value is found by the index of the gradient table //Here, the operands between Perm have no effect, just the color of the result c[di][dj][dk] = ranfloat[perm_x[(i + di) & 255] & perm_y[(j + dj) & 255] & perm_z[(k + dk) & 255]]; } } } return trilinear_interp(c, u, v, w); }};
- 输出:
从上图可以看出还是有明显的格子,其中有一部分是
mach band
效应,这是由线性变化的颜色构成的视觉感知效果。解决这一问题我们用:hermite cube
来平滑差值什么是mach bands效应?
简单来说,mach band效应是一种
边缘对比
效应。当观察两块亮度不同的区域时,边界亮度对比加强,使得轮廓更为明显mach bands效应的原因
人类的视觉系统可以增强边缘对比度
什么是hermite cube?
就是我们在前面提到柏林噪声时说的替换t,\(x = x1 + t(x2 – x1),t = 3t^2 – 2t^3\)
double noise(const vec3& p) const{ auto u = p.x() - floor(p.x()); auto v = p.y() - floor(p.y()); auto w = p.z() - floor(p.z()); u = 3 * std::pow(u, 2) - 2 * std::pow(u, 3); v = 3 * std::pow(v, 2) - 2 * std::pow(v, 3); w = 3 * std::pow(w, 2) - 2 * std::pow(w, 3); //...}
输出:
加入花纹
加入一个变量使得点的采样频率更高,实际上是
索引值&255也有周期
class noise_texture : public texture{ private: perlin noise; double scale; public: noise_texture() = default; noise_texture( double sc ) : scale(sc) {} virtual vec3 value(double u, double v, const vec3& p) const override { return vec3(1, 1, 1) * noise.noise(p * scale); }};hittable_list two_perlin_spheres(){ auto pertext = make_shared(10.0);}
- 输出:
上图现在看起来还是有格子,这是因为实现的max和min总是精确地到整数。因此我们将使用
梯度向量
,并以点乘的方法将min和max的值退离网格点(也就是变为小数)将参考点到小数部分(u,v,w)的向量和伪随机单位向量进行点乘,再进行线性插值
//vec3的三线性插值inline double vec3_trilinear_interp(vec3 c[2][2][2], double u, double v, double w){ auto uu = 3 * std::pow(u, 2) - 2 * std::pow(u, 3); auto vv = 3 * std::pow(v, 2) - 2 * std::pow(v, 3); auto ww = 3 * std::pow(w, 2) - 2 * std::pow(w, 3); auto lerp_x_y0_z0 = dot(vec3(u, v, w), c[0][0][0]) + (dot(vec3(u - 1, v, w), c[1][0][0]) - dot(vec3(u, v, w), c[0][0][0])) * uu; auto lerp_x_y1_z0 = dot(vec3(u, v - 1, w), c[0][1][0]) + (dot(vec3(u - 1, v - 1, w), c[1][1][0]) - dot(vec3(u, v - 1, w), c[0][1][0])) * uu; auto lerp_x_y0_z1 = dot(vec3(u, v, w - 1), c[0][0][1]) + (dot(vec3(u - 1, v, w - 1), c[1][0][1]) - dot(vec3(u, v, w - 1), c[0][0][1])) * uu; auto lerp_x_y1_z1 = dot(vec3(u, v - 1, w - 1), c[0][1][1]) + (dot(vec3(u - 1, v - 1, w - 1), c[1][1][1]) - dot(vec3(u, v - 1, w - 1), c[0][1][1])) * uu; auto lerp_y_z0 = lerp_x_y0_z0 + (lerp_x_y1_z0 - lerp_x_y0_z0) * vv; auto lerp_y_z1 = lerp_x_y0_z1 + (lerp_x_y1_z1 - lerp_x_y0_z1) * vv; auto lerp_z = lerp_y_z0 + (lerp_y_z1 - lerp_y_z0) * ww; return lerp_z;}class perlin{ private: vec3* ran_vec3;//Get a random vec3 //... public: perlin() { ran_vec3 = new vec3[point_count]; for (int i = 0; i < point_count; ++i) { ran_vec3[i] = unit_vector(vec3::random(-1,1)); } //... } ~perlin() { delete[] ran_vec3; //... } //Calculate the random value of position P double noise(const vec3& p) const { //decimal part auto u = p.x() - floor(p.x()); auto v = p.y() - floor(p.y()); auto w = p.z() - floor(p.z()); //integral part int i = floor(p.x()); int j = floor(p.y()); int k = floor(p.z()); vec3 c[2][2][2]; //xor //Record the Hasche transform to get the gradient of the reference point //The limiting gradient value is [0,255] for (int di = 0; di < 2; ++di) { for (int dj = 0; dj < 2; ++dj) { for (int dk = 0; dk < 2; ++dk) { //hash transform //Finds the index of the gradient table by the index values of the permutation table //The gradient value is found by the index of the gradient table //Here, the operands between Perm have no effect, just the color of the result c[di][dj][dk] = ran_vec3[perm_x[(i + di) & 255] ^ perm_y[(j + dj) & 255] ^ perm_z[(k + dk) & 255]]; } } } return vec3_trilinear_interp(c, u, v, w); }};//The output of the Berlin interpolation may be negative,These negative numbers become NaN when Gama corrects them by sqrt()//We should map the output to between 0 and 1class noise_texture : public texture{ //... public: virtual vec3 value(double u, double v, const vec3& p) const override { return vec3(1, 1, 1) * ( 1.0 + noise.noise( p * scale) ) * 0.5; }};
- 输出:
噪声叠加(扰动,turblence)
使用多个频率相加得到复合噪声
关于倍频系数weight的算法:这个weight是会不断变化的(否则并不会那么自然,也就是起伏过于平缓)
\[\LARGE frequency = 2^n\\\LARGE amplitude = persistence^i\\其中,frequency为频率,amplitude为振幅\]
class perlin{ public: //turbulence double turb(const vec3& p, int depth = 7)const { auto accum = 0.0; vec3 temp_p = p; auto weight = 1.0; for (int i = 0; i < depth; ++i) { accum += weight * noise(temp_p); weight *= 0.5;//amplitude temp_p *= 2;//frequence } return std::fabs(accum); }}
- 输出:
大理石纹理
调整颜色:让颜色与sin函数的值成比例
条纹起伏波荡:并使用turblence函数去调整相位(平移sin函数的参数),使得带状条纹起伏波荡
可以看到上图的频率很高,若我们只是降低频率是否可以实现大理石纹理呢?
class noise_texture : public texture{ //... public: virtual vec3 value(double u, double v, const vec3& p) const override { return vec3(1, 1, 1) * (1 + 5 * noise.turb(p)) * 0.5; }};
- 输出:可以看到并没有大理石纹理,亮度特别高,颜色值不够平缓
加入sin函数
class noise_texture : public texture{ //... public: virtual vec3 value(double u, double v, const vec3& p) const override { return vec3(1, 1, 1) * ( sin( 5 * noise.turb(p) ) ) * 0.5; }};
- 输出:可以看到还差点亮度,应该更白
实现大理石纹理
class noise_texture : public texture{ //... public: virtual vec3 value(double u, double v, const vec3& p) const override { return vec3(1, 1, 1) * ( 1 + sin( scale*p.z() + 10 * noise.turb(p) ) ) * 0.5; }};
纹理贴图
在这之前我们实现了程序式纹理,我们也能使用图像来进行改变物体外观。也就是,读取一张图片,将2D(u,v)坐标系映射在图片上
纹理贴图可以通过一个通常的纹理管线来描述:
如下对上图进行概述:
- 通过
投影方程
将空间
中的点转换为一组参数空间值
(也就是u,v坐标),这个值和纹理有关 - 使用
映射函数
将参数空间值转换到纹理空间
- 使用纹理空间值从
纹理
中获取相应的值
- 使用
值变换函数
对检索结果进行值变换,再使用得到的新值改变表面属性
- 来看个例子:
- 找到物体在空间中的位置(x,y,z)
- 使用投影函数转换为二元向量(u,v)
- 我们假设图像分辨率为(256*256),则分别乘以(u,v)上两个维度的坐标,得到纹理坐标
- 通过纹理坐标,在纹理贴图上查找到对应颜色值
- 假设找到的颜色值太暗,我们可以使用值变换函数,对每个向量乘以1.1
- 最后,我们就可以将作用于值变换函数的向量 用于着色方程
- 通过
映射函数的实现:
使用(u,v)的直接想法:
将u和v调整比例后取整,再将其对应到像素坐标(i,j)上
可行,但比较麻烦。试想每次图片分辨率发生变化,我们都要修改代码
更广泛的做法:
采用纹理坐标系代替图像坐标系,范围是[0,1]
如:对一般的图像来说,宽\(N_x\),高\(N_y\),他的一像素(i,j),在纹理坐标系中的坐标为:$$\large u = \frac{i}{N_x – 1}\
\large v = \frac{j}{N_y – 1}$$
球体uv坐标:
对于球体来说,uv的计算是基于球面坐标的,只需按比例转化即可得到uv坐标
当我们有一个球面坐标:\(\large (\theta,\phi)\),其中\theta是朝下距离极轴的角度,\phi是绕极轴旋转的角度将球面坐标映射到[0,1]:\(\large u = \frac{\phi}{2 \pi} v = \frac{\theta}{\pi}\)
计算\(\large \theta 和 \phi\)
将球面坐标系转化为直角坐标系:
\[\large x = cos(\theta) cos(\theta) \\\large y = sin(\theta) cos(\phi) \\\large z = sin(\theta)\]
将第二步倒过来,根据将直角坐标系转化为球面坐标系
利用库函数
atan2()
,给出任意一角度的sin和cos值,即可得出这个角的角度值//inverse functionstheta = atan2(y,x);//range [-pi, pi]phi = atan2(z);//range [-pi/2, pi/2]
实现工具函数
void get_sphere_uv(const vec3& p, double& u, double& v){ auto phi = std::atan2(p.z(), p.x()); auto theta = std::asin(p.y()); u = 1 - (phi + pi) / (2 * pi);//Consider phi + pi = 2pi v = (theta + pi / 2) / pi;}
更新sphere::hit():
//spere::hitget_sphere_uv( (rec.p - center) / radius, rec.u, rec.v );
存放图片
我们使用
stb_image.h
——将图片信息读入一个unsigned char(8bit)的数组中//Texture mappingclass image_texture : public texture{ private: unsigned char* data; int nx, ny; public: image_texture() = default; image_texture(unsigned char* pixels, int A, int B) : data(pixels), nx(A), ny(B) {} ~image_texture() { delete data; } virtual vec3 value(double u, double v, const vec3& p) const override { if (data == nullptr) { return vec3(0, 1, 1); } //texture-space locations auto i = static_cast(u * nx); auto j = static_cast((1 - v) * ny - 0.001); //because v doesn't take into account -2/Pi if (i < 0) i = 0; if (j nx - 1) i = nx - 1; if (j > ny - 1) j = ny - 1; auto r = static_cast(data[3 * i + 3 * nx * j + 0] / 255.0); auto g = static_cast(data[3 * i + 3 * nx * j + 1] / 255.0); auto b = static_cast(data[3 * i + 3 * nx * j + 2] / 255.0); return vec3(r, g, b); }};//main.cpp#define STB_IMAGE_IMPLEMENTATION#include "stb_image.h"hittable_list texture_mapping(){int nx, ny, nn;unsigned char* texture_data = stbi_load("Bronya.jpg", &nx, &ny, &nn, 0);auto bronya_surface = make_shared(make_shared(texture_data, nx, ny));auto bronya = make_shared(vec3(0, 0, 0), 3, bronya_surface);return hittable_list(bronya);}
输出:
这是我的素材这是我的输出结果
矩阵和光源发射光线的材质
- 加入发射函数。这个材质只要考虑自己发射的光线的颜色,无需考虑反射折射
class material{ public: //scatter virtual bool scatter(const ray& r_in, const hit_record& rec, vec3& attenuation, ray& scattered) const = 0; //emit virtual vec3 emitted(double u, double v, const vec3& p) const { return vec3(0, 0, 0); }};class diffus_light : public material{ private: shared_ptr emit; public: diffus_light(shared_ptr a) : emit(a) {} virtual bool scatter(const ray& r_in, const hit_record& rec, vec3& attenuation, ray& scattered) const override { return false; } virtual vec3 emitted(double u, double v, const vec3& p) const { return emit->value(u, v, p); }};
实现一个纯黑背景,并让所有光线都来自光源材质
在ray_color()中加入背景色变量,并由emitted()递归产生新颜色值
vec3 ray_color(const ray& r, const vec3& background, const hittable& world, int depth){ hit_record rec; if (depth emitted(rec.u, rec.v, rec.p); if (rec.mat_ptr->scatter(r, rec, attenuation, scattered) == false) { return emitted; } return emitted + attenuation * ray_color(scattered, background, world, depth - 1);}void output_image(){ //... const vec3 background(0, 0, 0); //... color += ray_color(r, background, world, max_depth);}
加入矩形光源
实现轴对齐的矩形
将矩形放在xy平面
使用z值确定平面,xy确定平面范围
。如:若z = k,一个轴对齐的矩形由\(\large x = x_0, x = x_1, y = y_0, y = y_1\)四条直线构成判断光线是否于矩形相交
\[\large 射线的z值为\; z(t) = a_z + t · \overrightarrow{b} \\\large 而z = k时 t的值为\; t = \frac{k – a_z}{\overrightarrow{b_z}} \\\large 即可求解x和y:\; \\\large x = a_x + t · \overrightarrow{b_x} \\\large y = a_y + t · \overrightarrow{b_y} \\\large 其中,若x在[x_0, x_1]且y在[y0,y1],射线就击中了矩阵\]
//xyclass xy_rect : public hittable{ private: double x0, x1, y0, y1, k; shared_ptr mp; public: xy_rect() = default; xy_rect( double _x0, double _x1, double _y0, double _y1, double _k, shared_ptr mat) : x0(_x0), x1(_x1), y0(_y0), y1(_y1), k(_k), mp(mat) {} virtual bool hit(const ray& r, const double t_min, const double t_max, hit_record& rec) const override; virtual bool bounding_box(double t0, double t1, aabb& output_box) const override;};bool xy_rect::hit(const ray& r, const double t_min, const double t_max, hit_record& rec) const{ auto t = (k - r.get_origin().z()) / r.get_direction().z(); if (t t_max) { return false; } auto x = r.get_origin().x() + t * r.get_direction().x(); auto y = r.get_origin().y() + t * r.get_direction().y(); if (x x1 || y y1) { return false; } rec.u = (x - x0) / (x1 - x0); rec.v = (y - y0) / (y1 - y0); rec.t = t; vec3 outward_normal = vec3(0, 0, 1); rec.set_face_normal(r, outward_normal); rec.mat_ptr = mp; rec.p = r.at(t); return true;}bool xy_rect::bounding_box(double t0, double t1, aabb& output_box) const{ output_box = aabb(vec3(x0, y0, k - 0.0001), vec3(x1, y1, k + 0.0001)); return true;}//main.cpphittable_list simple_light(){ hittable_list objects; auto pertext = make_shared(4); objects.add(make_shared(vec3(0, -1000, 0), 1000, make_shared(pertext))); objects.add(make_shared(vec3(0, 2, 0), 2, make_shared(pertext))); auto diff_light = make_shared(make_shared(vec3(4, 4, 4))); objects.add(make_shared(vec3(0, 7, 0), 2, diff_light)); objects.add(make_shared(3, 5, 1, 3, -2, diff_light)); return static_cast(make_shared(objects, 0.0, 1.0));}
输出:
- 球形光源:
- 矩形光源:
Cornell Box
xz平面和yz平面
//xzclass xz_rect : public hittable{ private: double x0, x1, z0, z1, k; shared_ptrmp; public: xz_rect() = default; xz_rect(double _x0, double _x1, double _z0, double _z1, double _k, shared_ptr mat) : x0(_x0), x1(_x1), z0(_z0), z1(_z1), k(_k), mp(mat) {} virtual bool hit(const ray& r, const double t_min, const double t_max, hit_record& rec) const override; virtual bool bounding_box(double t0, double t1, aabb& output_box) const override;};bool xz_rect::hit(const ray& r, const double t_min, const double t_max, hit_record& rec) const{ auto t = (k - r.get_origin().y()) / r.get_direction().y(); if (t t_max) { return false; } auto x = r.get_origin().x() + t * r.get_direction().x(); auto z = r.get_origin().z() + t * r.get_direction().z(); if (x x1 || z z1) { return false; } rec.u = (x - x0) / (x1 - x0); rec.v = (z - z0) / (z1 - z0); rec.t = t; vec3 outward_normal = vec3(0, 1, 0); rec.set_face_normal(r, outward_normal); rec.mat_ptr = mp; rec.p = r.at(t); return true;}bool xz_rect::bounding_box(double t0, double t1, aabb& output_box) const{ output_box = aabb(vec3(x0, k - 0.0001, z0), vec3(x1, k + 0.0001, z1)); return true;}//yzclass yz_rect : public hittable{ private: double y0, y1, z0, z1, k; shared_ptrmp; public: yz_rect() = default; yz_rect(double _y0, double _y1, double _z0, double _z1, double _k, shared_ptr mat) : y0(_y0), y1(_y1), z0(_z0), z1(_z1), k(_k), mp(mat) {} virtual bool hit(const ray& r, const double t_min, const double t_max, hit_record& rec) const override; virtual bool bounding_box(double t0, double t1, aabb& output_box) const override;};bool yz_rect::hit(const ray& r, const double t_min, const double t_max, hit_record& rec) const{ auto t = (k - r.get_origin().x()) / r.get_direction().x(); if (t t_max) { return false; } auto y = r.get_origin().y() + t * r.get_direction().y(); auto z = r.get_origin().z() + t * r.get_direction().z(); if (y y1 || z z1) { return false; } rec.u = (y - y0) / (y1 - y0); rec.v = (z - z0) / (z1 - z0); rec.t = t; vec3 outward_normal = vec3(1, 0, 0); rec.set_face_normal(r, outward_normal); rec.mat_ptr = mp; rec.p = r.at(t); return true;}bool yz_rect::bounding_box(double t0, double t1, aabb& output_box) const{ output_box = aabb(vec3(k - 0.0001, y0, z0), vec3(k + 0.0001, y1, z1)); return true;}
实现墙壁
hittable_list cornell_box(){ hittable_list objects; auto red = make_shared(make_shared(vec3(0.65, 0.05, 0.05))); auto white = make_shared(make_shared(vec3(0.73, 0.73, 0.73))); auto green = make_shared(make_shared(vec3(0.12, 0.45, 0.15))); auto light = make_shared(make_shared(vec3(15, 15, 15))); objects.add(make_shared(0, 555, 0, 555, 555, green)); objects.add(make_shared(0, 555, 0, 555, 0, red)); objects.add(make_shared(213, 343, 227, 332, 554, light)); objects.add(make_shared(0, 555, 0, 555, 0, white)); objects.add(make_shared(0, 555, 0, 555, 555, white)); objects.add(make_shared(0, 555, 0, 555, 555, white)); return static_cast(make_shared(objects, 0.0, 1.0));}void output_image(){ //... auto world = cornell_box(); vec3 lookfrom(278, 278, -800); vec3 lookat(278, 278, 0); vec3 up_vector(0, 1, 0); auto dist_to_focus = 10.0; auto aperture = 0.0; const auto aspect_ratio = double(image_width) / image_height; camera cam(lookfrom, lookat, up_vector, 40.0, aspect_ratio, aperture, dist_to_focus, 0.0, 1.0); //...}
输出:
加入两个长方体
Cornell Box里一般含有两个相对墙面来说有一定角度的长方体
我们先实现
轴对齐的长方体
,向其加入六个面
//Axis-aligned cuboidclass box : public hittable{ private: vec3 box_min; vec3 box_max; hittable_list sides; public: box() = default; box(const vec3& p0, const vec3& p1, shared_ptr ptr); virtual bool hit(const ray& r, const double t_min, const double t_max, hit_record& rec) const override { return sides.hit(r,t_min, t_max, rec); } virtual bool bounding_box(double t0, double t1, aabb& output_box) const override { output_box = aabb(box_min, box_max); return true; }};inline box::box(const vec3& p0, const vec3& p1, shared_ptr ptr){ box_min = p0; box_max = p1; sides.add(make_shared(p0.x(), p1.x(), p0.y(), p1.y(), p1.z(), ptr)); sides.add(make_shared(p0.x(), p1.x(), p0.y(), p1.y(), p0.z(), ptr)); sides.add(make_shared(p0.x(), p1.x(), p0.z(), p1.z(), p1.y(), ptr)); sides.add(make_shared(p0.x(), p1.x(), p0.z(), p1.z(), p0.y(), ptr)); sides.add(make_shared(p0.y(), p1.y(), p0.z(), p1.z(), p1.x(), ptr)); sides.add(make_shared(p0.y(), p1.y(), p0.z(), p1.z(), p0.x(), ptr));}
- 输出:
旋转&平移长方体
在光追中,我们通常用
instance
来完成旋转、平移
instance是一种经过旋转或平移等操作的几何图元
平移
我们并不需要去移动图元,我们只需对ray进行操作。如下图中将粉色盒子平移至右边的盒子,将位于原点的粉色盒子所有的组成部分的x+2。或者在hit函数中,把射线的原点-2计算得到的是左边的粉色,最后再+2得到右边盒子
class translate : public hittable{ private: shared_ptrptr; vec3 offset; public: translate(shared_ptr p, const vec3& displacement) : ptr(p), offset(displacement) {} virtual bool hit(const ray& r, const double t_min, const double t_max, hit_record& rec) const override; virtual bool bounding_box(double t0, double t1, aabb& output_box) const override;};bool translate::hit(const ray& r, const double t_min, const double t_max, hit_record& rec) const{ ray moved_r(r.get_origin() - offset, r.get_direction(), r.get_time()); if (ptr->hit(moved_r, t_min, t_max, rec) == false) { return false; } rec.p += offset; rec.set_face_normal(moved_r, rec.normal); return true;}bool translate::bounding_box(double t0, double t1, aabb& output_box) const{ if (ptr->bounding_box(t0, t1, output_box) == false) { return false; } output_box = aabb(output_box.get_min() + offset, output_box.get_max() + offset); return true;}
旋转
从零手写一个软渲染器day3 MVP变换 – 爱莉希雅 – 博客园 (cnblogs.com)一些推导过程
绕z轴旋转
\[\large x’ = cos(\theta)·x – sin(\theta)·y \\\large y’ = sin(\theta)·x + cos(\theta)·y\]
绕y轴旋转
\[\large x’ = cos(\theta)·x + sin(\theta) · z \\\large z’ = -sin(\theta) · x + cos(\theta)·z\]
绕x轴旋转
\[\large y’ = cos(\theta) · y – sin(\theta) · z \\\large z’ = sin(\theta) · y + cos(\theta) · z\]
旋转时,图元表面的法向量也发生了改变,因此还需重新计算法向量
class rotate_y : public hittable{ private: shared_ptr ptr; double sin_theta; double cos_theta; bool hasbox; aabb bbox; public: rotate_y(shared_ptr p, double angle); virtual bool hit(const ray& r, const double t_min, const double t_max, hit_record& rec) const override; virtual bool bounding_box(double t0, double t1, aabb& output_box) const override { output_box = bbox; return hasbox; }};//Calculate the coordinates of the six rotating points of the axially aligned cubeinline rotate_y::rotate_y(shared_ptr p, double angle) : ptr(p){ auto radians = degrees_to_radians(angle); sin_theta = sin(radians); cos_theta = cos(radians); hasbox = ptr->bounding_box(0, 1, bbox); vec3 min(infinity, infinity, infinity); vec3 max(-infinity, -infinity, -infinity); for (int i = 0; i < 2; ++i) { for (int j = 0; j < 2; ++j) { for (int k = 0; k < 2; ++k) { auto x = i * bbox.get_max().x() + (1 - i) * bbox.get_min().x(); auto y = j * bbox.get_max().y() + (1 - j) * bbox.get_min().y(); auto z = k * bbox.get_max().z() + (1 - k) * bbox.get_min().z(); auto newx = cos_theta * x + sin_theta * z; auto newz = -sin_theta * x + cos_theta * z; vec3 tester(newx, y, newz); for (int c = 0; c hit(rotated_r, t_min, t_max, rec) == false) { return false; } vec3 normal = rec.normal; vec3 p = rec.p; p[0] = cos_theta * rec.p[0] + sin_theta * rec.p[2]; p[1] = -sin_theta * rec.p[0] + cos_theta * rec.p[2]; normal[0] = cos_theta * rec.normal[0] + sin_theta * rec.normal[2]; normal[0] = -sin_theta * rec.normal[0] + cos_theta * rec.normal[2]; rec.p = p; rec.set_face_normal(rotated_r, normal); return true;}
输出:
Volume
volume rendering(体渲染):用于显示离散三维采样数据集的二维投影的技术。主要用于渲染
体积云、体积光
将
volume表示为一个随机表面
加入volume部分内容会导致代码结构混乱,因为volume和平面表面是完全不同的两种东西。但我们可以将volume表示为一个随机表面。如一团烟雾可以用在概率上不确定位置的平面来代替实现固定密度的volume
- ray可以在volume内部任意一点发生散射,也可以直接穿过volume,volume越薄,直接穿过的概率更大;volume越浓,越可能发生散射
- 在任意微小的距离差\(\large \Delta L\)发生散射的概率:
- \[\large probability = C· \Delta L \\其中,C是volume的光学密度比例常数\]
- 经过一系列不同的等式运算,将会随机得到一个光线发生散射的距离值。根据这个距离值,若散射点在volume外部,则认为没有相交
- 对于静态的体积体,我们只需它的密度C和边界
//The boundary of a volumeclass constant_medium : public hittable{ private: shared_ptr boundary; shared_ptr phase_function; double neg_inv_density; constant_medium() = default; public: constant_medium(shared_ptr b, double d, shared_ptr a) : boundary(b), neg_inv_density(-1 / d) { phase_function = make_shared(a); } virtual bool hit(const ray& r, const double t_min, const double t_max, hit_record& rec) const override; virtual bool bounding_box(double t0, double t1, aabb& output_box) const override { return boundary->bounding_box(t0, t1, output_box); }};bool constant_medium::hit(const ray& r, const double t_min, const double t_max, hit_record& rec) const{ //Print samples when debugging. To enable, set enableDebug true. const bool enableDebug = false; const bool debugging = enableDebug && random_double() hit(r, -infinity, infinity, rec1) == false) { return false; } if (boundary->hit(r, rec1.t + 0.0001, infinity, rec2) == false) { return false; } if (debugging == true) { std::cerr << std::endl << "t0 = " << rec1.t << ", t1 = " << rec2.t << std::endl; } if (rec1.t t_max) rec2.t = t_max; if (rec1.t >= rec2.t) return false; if (rec1.t distance_inside_boundary) { return false; } rec.t = rec1.t + hit_distance / ray_length; rec.p = r.at(rec.t); if (debugging) { std::cerr << "hit_distance = " << hit_distance << '\n' << "rec.t = " << rec.t << '\n' << "rec.p = " << rec.p << '\n'; } rec.normal = vec3(1, 0, 0); rec.front_face = true; rec.mat_ptr = phase_function; return true;}class isotropic : public material{ private: shared_ptr albedo; public: isotropic(shared_ptr a) : albedo(a) {} virtual bool scatter(const ray& r_in, const hit_record& rec, vec3& attenuation, ray& scattered) const override { scattered = ray(rec.p, random_in_unit_sphere(), r_in.get_time()); attenuation = albedo->value(rec.u, rec.v, rec.p); return true; }};
最终成品
- 将所有实现的内容渲染出来,并用体积雾盖住,再加入一个在电解质内部充满volume的蓝色球体
static hittable_list final_scene() { hittable_list boxes1; auto ground = make_shared(make_shared(vec3(0.48, 0.83, 0.53))); const int boxes_per_side = 20; for (int i = 0; i < boxes_per_side; i++) { for (int j = 0; j < boxes_per_side; j++) { auto w = 100.0; auto x0 = -1000.0 + i * w; auto z0 = -1000.0 + j * w; auto y0 = 0.0; auto x1 = x0 + w; auto y1 = random_double(1, 101); auto z1 = z0 + w; boxes1.add(make_shared(vec3(x0, y0, z0), vec3(x1, y1, z1), ground)); } } hittable_list objects; objects.add(make_shared(boxes1, 0, 1)); auto light = make_shared(make_shared(vec3(7, 7, 7))); objects.add(make_shared(123, 423, 147, 412, 554, light)); auto center1 = vec3(400, 400, 200); auto center2 = center1 + vec3(30, 0, 0); auto moving_sphere_material = make_shared(make_shared(vec3(0.7, 0.3, 0.1))); objects.add(make_shared(center1, center2, 0, 1, 50, moving_sphere_material)); objects.add(make_shared(vec3(260, 150, 45), 50, make_shared(1.5))); objects.add(make_shared( vec3(0, 150, 145), 50, make_shared(vec3(0.8, 0.8, 0.9), 10.0) )); auto boundary = make_shared(vec3(360, 150, 145), 70, make_shared(1.5)); objects.add(boundary); objects.add(make_shared( boundary, 0.2, make_shared(vec3(0.2, 0.4, 0.9)) )); boundary = make_shared(vec3(0, 0, 0), 5000, make_shared(1.5)); objects.add(make_shared( boundary, .0001, make_shared(vec3(1, 1, 1)))); int nx, ny, nn; auto tex_data = stbi_load("Bronya.jpg", &nx, &ny, &nn, 0); auto emat = make_shared(make_shared(tex_data, nx, ny)); objects.add(make_shared(vec3(400, 200, 400), 100, emat)); auto pertext = make_shared(0.1); objects.add(make_shared(vec3(220, 280, 300), 80, make_shared(pertext))); hittable_list boxes2; auto white = make_shared(make_shared(vec3(0.73, 0.73, 0.73))); int ns = 1000; for (int j = 0; j < ns; j++) { boxes2.add(make_shared(vec3::random(0, 165), 10, white)); } objects.add(make_shared( make_shared( make_shared(boxes2, 0.0, 1.0), 15), vec3(-100, 270, 395) ) ); return static_cast(make_shared(objects, 0.0, 1.0));}
- 输出:借用一下原教程的图,因为太慢了!
reference
Ray Tracing: The Next Week V3.0中文翻译(上) – 知乎 (zhihu.com)
QianMo/Real-Time-Rendering-3rd-CN-Summary-Ebook: 电子书 -《Real-Time Rendering 3rd》提炼总结 | 全书共9万7千余字。你可以把它看做中文通俗版的《Real-Time Rendering 3rd》,也可以把它看做《Real-Time Rendering 3rd》的解读版与配套学习伴侣,或者《Real-Time Rendering 4th》的前置阅读材料。 (github.com)
GAMES101:现代计算机图形学入门 – 计算机图形学与混合现实在线平台 (games-cn.org)
Wikipedia, the free encyclopedia (jinzhao.wiki)
[Nature of Code] 柏林噪声 – 知乎 (zhihu.com)
(26条消息) 什么是马赫带效应_李守聪的博客-CSDN博客_马赫带效应
3D Math Primer for Graphics and Game Development 2nd
gpu gems1