inksoul
/
InkSoul
1
0
Fork 0
Code Issues Pull Requests Packages Projects Releases Wiki Activity
InkSoul/content/computergraphic/路径追踪.md

537 lines
17 KiB
Markdown
Raw Permalink Normal View History

init
2023-06-18 23:16:41 +08:00
---
title: "路径追踪"
date: 2022-07-24T15:01:07+08:00
---
### 概率论基础复习
在路径追踪的计算中需要引入部分概率论内容，在此进行部分复习
#### 随机变量
* 随机变量：一个值，表示某些事件发生
* $x \backsim p(x)$ ：随机变量的分布，随机变量取到某个值的可能性
例如：对于一个骰子，共有6个面，每个面在上面的概率相等，都为1/6
#### 概率
* 概率必定为非负
* 所有事件的发生概率和为1
$$p_i \geq 0$$
$$\sum_{i = 1}^n p_i = 1$$
例如：对于上述的骰子，概率为$p_i = \frac{1}{6}$
#### 随机变量的期望
如图，对于这四个可供选择的长方形，每一个选择有一个值为$X_i$,每一个选择都对应一个概率为$p_i$,那么期望则为将所有值乘以其对应的概率后累加得到的数值
![](../../images/eg_expected_value.png)
可见上图的期望值为
$$E[X] = \sum_{i = 1} ^ n x_ip_i$$
例如：一个骰子可以取值1、2、3、4、5、6,每一个取值的概率都是1/6,期望为
$$E[X] = \sum_{i=1}^n \frac{i}{6} = (1+2+3+4+5+6)/6 = 3.5$$
#### 概率密度函数
随机变量的取值可以是一些固定的值(离散型)，也可以是下图所示一些连续的值(连续型)
![](../../images/probability_distribution_function.png)
某一位置的概率为其周围一小段微元与曲线相连线形成的梯形面积
上面这条曲线P(x)即为概率密度函数(简称PDF)
概率密度在所有值上的积分等于1，与概率和为1类似
期望的计算则是用任意一个取值乘以概率密度，并加以积分
对于满足下列条件的P(x)
$$p(x) \geq 0 \\ \\ and \int p(x)dx = 1$$
期望值为
$$E[X] = \int x P(x)dx$$
当一个函数本身是随机变量时：
例如一个随机变量满足一定的PDF
$$X \backsim p(x)$$
另外一个函数为
$$Y = f(X)$$
此时Y的期望为
$$E[Y] = E[f(X)] = \int f(x)p(x)dx$$
--------------------
### 蒙特卡洛积分(Monte Carlo Integration)
#### 使用理由
计算给定函数的定积分时，如果$f(x) = x^2$,就可以计算出不定积分为$\frac{1}{3}X^3 + 常数$，之后利用不定积分即可求出它在a和b的值，两者相减即得积分
若函数较为复杂，如下图所示，不便于使用解析方法，就要使用数值的方法计算(蒙特卡洛积分)
![](../../images/reason_for_monte_carlo_integration.png)
#### 大致流程
黎曼积分：
将上图a和b之间均匀拆分为100份，取每一份的中间位置，找到对应的Y(即将整个曲线下方面积分解为各个小长方形面积之和)
蒙特卡洛积分：
* 在a和b之间随便取一个数，找到这个数对应的f(x)
* 假设整个曲线为一个长方形，长方形高度为刚才去的值的对应位置高度，长方形宽度视为a到b的长度，利用长方形面积近似曲线下围面积
* 将第二步重复多次，在a到b区间采样多次，将得到的长方形面积平均求值，就可以得到相对准确的结果
#### 定义
定积分f(x)，即从a到b的值：
$$\int_a^b f(x)dx$$
x的积分域在ab区间，对概率密度函数随机采样一个变量，得到Xi
$$随机变量 \\ \\ \\ X_i \backsim p(x)$$
蒙特卡洛积分可近似为下列式子，即$f(Xi)/p(xi)$的平均
$$F_N = \frac{1}{N} \sum_{i = 1}^N {\frac{f(X_i)}{p(X_i)}}$$
#### 特殊情况
* 在a到b之间均匀采样时，我们使用的概率密度函数(PDF)在各处的值相等，即为一个常数，又已知PDF在积分域上的积分值为1，此PDF可解出值为$\frac{1}{(b-a)}$
![](../../images/uniform_monte_carlo_estimator.png)
* 采用蒙特卡洛积分计算时，需要知道f(Xi)和p(Xi),随机采样得到的值为Xi,蒙特卡洛积分只需计算f(Xi)除以p(Xi)的平均值，p(Xi)的值始终为$\frac{1}{(b-a)}$
$$F_N = \frac{b-a}{N}\sum_{i=1}^Nf(X_i)$$
#### 总结
蒙特卡洛积分：
只需积分域内使用一个PDF采样，得到该样本的f(Xi)和概率密度P(Xi)的值，两者相除后取平均值
注意：
1. N越大(采样次数越多)，得到的结果越精准
2. 积分域是定义在X上的积分，所以采样目标一定是X
### 路径追踪(path tracing)
#### Whitted-Style Ray Tracing 存在的问题
* Whitted-style Ray Tracing 对Specular(镜面)材质的表现效果是正确的，而glossy(磨砂)材质的效果则表现不正确
Specular材质：光线在物体表面能发生镜面反射则称为此材质，光线在打在表面时，会沿着镜面反射方向离去，也是镜面能够映出周围环境光的原因
Glossy材质：有镜面反射的样子，但又有些许模糊，类似毛玻璃效果，有一定粗糙度但能产生高光
![](../../images/mirror_and_glossy_reflection.png)
* whitted-style ray tracing 在遇到漫反射时会停止光线弹射直接做相应的shading
这种情况会忽视漫反射物体和漫反射物体之间的光线都不被渲染，也不符合漫反射将光线均匀地反射到各个不同方向上的物理定义
两个示例中左侧为直接光照，右侧为全局光照，都采用路径追踪的方式计算得到
![](../../images/path_traced_direct_global_illumination.png)
可见在场景中只有上面一个光源的情况下，在直接光照里，天花板为黑色，但实际上光线会发生多次弹射，天花板应是亮的，右侧才是我们想要的效果
Color bleeding: 左侧全局光照图中高立方体左侧是红，光线在弹射到墙面后再反射到这个面，之后到达摄像头(眼睛)，z这种现象就为Color bleeding
Cornell box真实存在，3D模型也是存在的，广泛用于测试各种全局光照效果
----------------------
两中表现效果错误来自于Whitted-style ray tracing的光线弹射计算
1. 光线打到specular的物体上，会沿着镜面方向反射或沿着折射方向折射
2. 光线打到diffuse(漫反射)的物体上，这条光线就停止了
虽然存在错误，但完全按照物理量推算的渲染方程是完全正确的，为了正确计算，我们要解出这个渲染方程
$$L_o(p,w_o) = L_e(p,w_o) + \int_{\Omega+} {L_i(p,w_i)f_r(p,w_i,w_o)(n\cdot w_i)dw_i}$$
#### 蒙特卡洛积分分解渲染方程
##### 直接光照
假设存在下图场景，我们认为各个方向进来的光$\omega i$ 均匀分布在球面上
![](../../images/suppose_scene_1.png)
* 假设该点不发光，因而直接光照强度结果来自于四面八方入射来的光照强度
忽略渲染方程发光项后的形式
$$L_o(p,w_o) =\int_{\Omega+} {L_i(p,w_i)f_r(p,w_i,w_o)(n\cdot w_i)dw_i}$$
###### 使用蒙特卡洛积分求解
* 通过蒙特卡洛求解我们需要在不同方向上采样，考虑对应的f(x)和PDF的值
着色点为p点时，从P点反射到摄像头的辐射为
$$L_o(p,w_o) =\int_{\Omega+} {L_i(p,w_i)f_r(p,w_i,w_o)(n\cdot w_i)dw_i}$$
蒙特卡洛为了算积分会在积分域上进行采样获得样本X，计算f(x)/p(x)，再求平均
$$\int_a^b f(x)dx \approx \frac{1}{N}\sum_{k=1}^N \frac{f(X_k)}{p(X_k)}$$
$$X_k \backsim p(x)$$
* 求出f(x)
将蒙特卡洛积分的方法迁移到解渲染方程上，则f(x)就是渲染方程中积分的所有计算式
$$L_i(p,w_i)f_r(p,w_i,w_o)(n\cdot w_i)$$
* PDF的值
PDF涉及对积分域进行采样，这里为对半球进行采样
使用最简单的均匀采样，认为半球上采样到任何一个方向上的概率密度是相同，PDF为一个常数
$$p(\omega_i) = \frac{1}{2\pi}$$
球面面积为$4\pi$，半球面面积为$2\pi$,半球对应立体角为$2\pi$,均匀采样使PDF的所在角度积分为1，因此PDF为$\frac{1}{2\pi}$
* 将渲染方程改写成蒙特卡洛形式
每次均匀地在半球上采样，取一个入射方向$\omega_i$,将上值代入可得下列式子
$$L_o(p,w_o) =\int_{\Omega+} {L_i(p,w_i)f_r(p,w_i,w_o)(n\cdot w_i)dw_i}$$
$$\approx \frac{1}{N}\sum_{i=1}^N \frac{L_i(p,w_i)f_r(p,w_i,w_o)(n\cdot w_i)}{p(\omega_i)}$$
由此可写出以下着色算法
```
着色点p，方向为w0
随机选择PDF采样得到的N个方向wi
初始化结果Lo
对于选中的wi
从p点连出一条光线
如果光线打到了光源
Lo+=(1/N) * Li * fr * cosine/pdf(wi)
返回Lo
```
-------------
```
shade(p, wo)
Randomly choose N directions wi~pdf
Lo = 0.0
For each wi
Trace a ray r(p, wi)
If ray r hit the light
Lo += (1 / N) * L_i * f_r * cosine / pdf(wi)
Return Lo
```
##### 全局光照
全局光照下，我们需要考虑反射面反射过来的光，计算从Q反射到P点反射了多少辐射能
![](../../images/introducing_global_illumination.png)
只需在直接光照的算法中添加分支便可支持全局光照
```
着色点p，方向为w0
随机选择PDF采样得到的N个方向wi
初始化结果Lo
对于选中的wi
从p点连出一条光线
如果光线打到了光源
Lo+=(1/N) * Li * fr * cosine/pdf(wi)
如果光线打到了q点的物体
Lo +=(1/N) * shade(q,-wi) * fr * cosine/pdf(wi)
返回Lo
```
--------------------------
```
shade(p, wo)
Randomly choose N directions wi~pdf
Lo = 0.0
For each wi
Trace a ray r(p, wi)
If ray r hit the light
Lo += (1 / N) * L_i * f_r * cosine / pdf(wi)
Else If ray r hit an object at q
Lo += (1 / N) * shade(q, -wi) * f_r * cosine / pdf(wi)
Return Lo
```
打到物体上是需要考虑q点反射过来的能量，即在q点-wi方向看过去的直接光照
---------------------
##### 上述算法存在的问题
* 这种方式打出不同光线再以递归的方式计算，会使光线的数量呈指数爆炸增长
如下图，打到第一个物体上后反射出N根光线，打到第二个物体后会再发出N个光线，即变成$N^2$条光线，会超出处理能力
![](../../images/explosion_of_ray_as_ray_go_up.png)
解决方法：
1. 取N等于1，但是会产生较大噪声
2. 采用N=1来做蒙特卡洛积分，即为路径追踪，(N!=1时为分布式光线追踪，会产生指数爆炸)
3. N=1的情况下会产生噪声，但需要得到的只是一个像素的Radiance，所以只需计算穿过一个像素的多个path再求平均即可得到目标Radiance
path:一条完整连接视点和光源的路径
算法：
```
ray_generation(camPos, pixel)
Uniformly choose N sample positions within the pixel
pixel_radiance = 0.0
For each sample in the pixel
Shoot a ray r(camPos, cam_to_sample)
If ray r hit the scene at p
pixel_radiance += 1 / N * shade(p, sample_to_cam)
Return pixel_radiance
```
* 递归算法永远不停
真实世界中光线本身弹射次数也是不停的，但计算机不能无限模拟，也不能提前限制弹射次数，这样会损失多次弹射的能量
解决方案：
用Russian Roulette(俄罗斯轮盘赌)的方法来决定以一定概率去停止光线继续弹射
俄罗斯轮盘思想：
1. 某一着色点出射的Radiance是$L_o$
2. 自定一个概率p(0< p < 1)
3. 以一定的概率p往某一方向打一条光线
4. 得到一定结果后取$L_o/p$为返回值
5. 在1-p的概率内就不打光线，返回值为0
可视为取两个值的离散型随机变量，可以通过概率乘以值并加起来计算期望
$$E = P*(L_o/p) + (1-p)*0 = L_o$$
对应修改后的着色算法
指定概率$P_RR$
随意在0到1中取一个数ksi
如果$ksi>P_RR$,意味着不该往外把一根光线
其他情况下正常打出光线，最后结果需除以$P_RR$
```
shade(p, wo)
Manually specify a probability P_RR
Randomly select ksi in a uniform dist. in [0, 1]
If (ksi > P_RR) return 0.0;
Randomly choose ONE direction wi~pdf(w)
Trace a ray r(p, wi)
If ray r hit the light
Return L_i * f_r * cosine / pdf(wi) / P_RR
Else If ray r hit an object at q
Return shade(q, -wi) * f_r * cosine / pdf(wi) / P_RR
```
至此得出不太高效的正确path tracing算法
----------------
Samples(采样率)可认为是path，SPP(samples per pixel)就是一个像素打出多少path
low SPP下，计算速度快，图片噪声多
High SPP下，计算速度慢，图片噪声少
![](../../images/High_and_low_spp_result.png)
##### 算法依旧不高效的原因
![](../../images/reason_of_being_inefficient.png)
光线能否打到光源上取决于运气，当光源小时往往会浪费过多光线，需要改用更好的PDF用来采样，而非简单的均匀采样
##### 光源上采样
蒙特卡洛允许许多采样方法，我们可以直接在光源上采样
![](../../images/sample_the_light.png)
n'为光源本身朝向，在与着色点连线后可得$\theta$(连线和着色点法线的夹角) $\theta$'(连线和光源法线的夹角)，此时若在视为二维平面的光源上均匀采样，则PDF为$\frac{1}{A}$
此时得到一个定义在光源上的积分，需要将渲染方程改写
* 将对$d\omega$的积分改为对$dA$的积分
dA: 光源上的一个小平面
$d\omega$：是上面小的表面投影到单位球上形成的面积(单位球半径为1，立体角为dA投影到单位球上的面积除以1的平方)
$d\omega$和dA的关系：dA的朝向不一定朝向着色点，则需要计算$dAcos\theta'$，将光源所在平面垂直地和着色点连线，除以距离绝对值的平方，得出$d\omega$
$$d\omega = \frac{dAcos\theta'}{||x' - x||^2}$$
重写后的渲染方程为：
$$L_o(x,\omega_o) = \int_{\Omega+} L_i(x,\omega_i)f_r(x,\omega_i,\omega_o)\cos \theta d\omega_i$$
$$\int_A L_i(x,\omega_i)f_r(x,\omega_i,\omega_o)\frac{dAcos\theta'}{||x' - x||^2} dA$$
即可使用蒙特卡洛积分进行计算,此时的PDF为1/A
##### 改进后的算法
着色点的着色结果来源于两个部分：
* 光源的贡献(直接光源，无需进行俄罗斯轮盘赌剔除)
均匀地在光源上采样后对改写后的渲染方程使用蒙特卡洛积分计算
* 非光源贡献(间接光源，需要进行俄罗斯轮盘赌剔除)
如果通过了俄罗斯轮盘赌测试，就发出一条射线，打到了一个非光源的点q时将其贡献加入
算法：
```
shade(p, wo)
# Contribution from the light source.
Uniformly sample the light at x’ (pdf_light = 1 / A)
L_dir = L_i * f_r * cos θ * cos θ’ / |x’ - p|^2 / pdf_light
# Contribution from other reflectors.
L_indir = 0.0
Test Russian Roulette with probability P_RR
Uniformly sample the hemisphere toward wi (pdf_hemi = 1 / 2pi)
Trace a ray r(p, wi)
If ray r hit a non-emitting object at q
L_indir = shade(q, -wi) * f_r * cos θ / pdf_hemi / P_RR
Return L_dir + L_indir
```
##### 处理光源遮挡
上述算法中，光源采样并未考虑光源遮挡问题
发生下图蓝色物体遮挡在光源和着色点之间时，需要判断光源能否贡献到着色点
![](../../images/light_source_blocked_judge.png)
判断方式：
着色点到光源上的采样点取一根连线，从着色点往这根连线的方向打出一根光线，判断是否会打到某个物体，自然可以直接判断直接光照是否被遮挡，不被遮挡时，直接光照可直接计算，被遮挡时，自然为0
改进后算法：
```
shade(p, wo)
# Contribution from the light source.
Uniformly sample the light at x' (pdf_light = 1 / A)
Shoot a ray from p to x’
If the ray is not blocked in the middle
L_dir = L_i * f_r * cos θ * cos θ’ / |x' - p|^2 / pdf_light
# Contribution from other reflectors.
L_indir = 0.0
Test Russian Roulette with probability P_RR
Uniformly sample the hemisphere toward wi (pdf_hemi = 1 / 2pi)
Trace a ray r(p, wi)
If ray r hit a non-emitting object at q
L_indir = shade(q, -wi) * f_r * cos θ / pdf_hemi / P_RR
Return L_dir + L_indir
```
最终path tracing的出的效果
![](../../images/is_path_tracing_correct.png)
可见，path tracing 与照片相比能够做到几乎百分百相似
##### 光线追踪新旧计算算法
早期多是Whitted-style Ray Tracing
现代：
可理解为所有光线传播方式的集合
* (Unidirectional & bidirectional) path tracing
* Photon mapping
* Metropolis light transport
* VCM / UPBP…
当然，现在还做不到完美的实时光线追踪，也是现在探讨前沿
----------------------------
##### 补充问题
* 蒙特卡洛积分应该选择什么样的PDF？
重要性采样理论
* 随机数是否有质量之分？
存在质量之分，不同的随机数有各自特点，如low discrepancy sequences、TAA中使用的Halton sequence
* 是否可以将采样半球和采样光源这两种采样方法结合起来，使其效果更好？
存在，称为multiple imp Sampling,简称MIS
* 我们算出了一个像素的Radiance，可是我们最后看到的是颜色，这个Radiance是不是颜色？
颜色与Radiance并非一一对应关系，将Radiance换算为颜色需要经过gamma校正