光子映射学习#

主要按[Realistic Image Synthesis Using Photon Mapping]这个书来学习光子映射

[A Practical Guide to Global Illumination using Photon Maps]大致是其精简版，也包括一些不同信息可以参考。

持续更新

传统光线追踪至少有两个问题。
- diffuse inter-reflection
没理解他说的是颜色不准还是指noise问题

caustics问题

这两个问题都是间接光问题。光子映射可以利用预计算的方法来解决。

是光线追踪的一个补充，能更有效地解决问题。

光子映射是一个two-pass方法。目的是计算diffuse面的间接光。
第一遍从光源发出光子，打到非镜面时存下信息。所有光子都在diffuse面上。
第二遍用统计方法，根据第一遍的光子映射信息算出radiance。

一个重点是光子映射与场景的geo信息是解耦的。所以可以用于非常复杂的场景。

1. Photon emission#

每个光子带相等的flux(辐射功率)。

点光源，光子的发射是均匀分布。
方相光源，从场景外均匀发射。
面光源？向半球发射，cos分布，normal方向概率大。

一个光源的总功率是 $P_{l i g h t}$ ，那么发射 $n_{e}$ 个光子，每个光子的功率

$P_{p h o t o n} = \frac{P_{l i g h t}}{n_{e}}$

多光源#

每个光源都要发射光子。

并不是多光源就需要更多光子，维持一个总量即可，各个光源按能量分配光子数。

还可使用importance sampling，对某些光源侧重。

2. Photon tracing#

photon tracing与ray tracing类似。
区别：光子追踪是散播flux，光线追踪是收集radiance。
因此与材质的交互也有不同。

一个例子是对于radiance，折射由折射率决定。光子不是这样算。光子打到物体上，可能反射，传播，或吸收。概率由material的反射率决定。

这里和ray tracing有巨大的区别。

我一开始没有仔细看这一节，而是直接用ray tracing的方法计算光子经过brdf处理后的能量，这样不断往下传播。也能出似是而非的图，但是有一系列问题。

论文实际需要用Russian roulette来决定对光子的操作。
一来可以有效的去除不重要的光子，减少计算量。
二来可以使每个光子的能量差不多。后续计算能量时效果更好。brdf处理后能量可能差距巨大。

反射还是吸收？#

设r为我们取的反射概率

p = randint() # [0-1]的随机数
if p < r:
    以原始能量反射
else:
    吸收

为什么这样操作是对的。比如材质r=0.7，有100个光子。

0.7的两种解读： 1. 反射100个光子，每个的能量为0.7。 2. 反射70个光子，每个的能量不变。

两种都是对的。我们就可以使用第二种解读。保持正确性的前提下，可以有效地减少计算。

反射，传播，还是吸收？#

以单色为例，设d为diffuse反射概率，s为specular反射概率。可以按三个概率区间决定光子如何传播。

[0, d] diffuse反射
[d, s + d] specular反射
[s + d, 1] 吸收

对于多色。

Photon storing#

只有光子碰到diffuse面时才储存。储存在photon map里。
每个光子在路径中可存多次。
存的数据为[位置，能量，入射角度]。

扩展到participating media#

使用光子映射能更容易地处理participating media。比如烟雾。

三个photon map#

为了效率，把光子分到三个map里。具体原因后续会讲。

Caustic photon map

在碰到diffuse面之前至少有过一次镜面反射的光子。 $L S^{+} D$
Global photon map

所有diffuse面的GI。 $L {S | D | V}^{*} D$
Volume photon map

participating media的间接光。 $L {S | D | V}^{+} V$

上面用了[Heckbert90]的描述光子路径的语法。

L 从光源发射
S 镜面反射或传输
D 非镜面反射或传输
V volume散射
{x|y|z} xyz中的一个
$x^{+}$ 1个或以上x
$x^{*}$ 任意个x

kd-tree#

需要频繁地查找一个光子附近的光子。需要适应非常不均匀分布的物体。

这里最常用的数据结构是平衡kd-tree。需要自己实现。

3. 计算面上的辐射率#

photon map可以看成输入功率的表示，计算radiance就是要对其积分。

点 $x$ 上的reflected radiance:

$L_{r} (x, ω) = \int_{Ω_{x}} f_{r} (x, ω^{'}, ω) L_{i} (x, ω^{'}) | c o s θ_{i} | d ω_{i}^{'}$

ω

为出射方向，

ω^{'}

为入射方向，

L_{i}

为入射radiance。

半球上所有方向(

ω^{'}

)的入射radiance，经过brdf和cos处理后积起来。

需要把radiance和flux的关系展开。

$E = \frac{d Φ}{d A}$

$L = \frac{d E}{c o s θ d ω} = \frac{d^{2} Φ}{c o s θ d A d ω}$

$L_{i} (x, ω^{'}) = \frac{d^{2} Φ_{i} (x, ω^{'})}{c o s θ_{i} d ω_{i}^{'} d A_{i}}$

代入 $L_{r}$

$L_{r} (x, ω) = \int_{Ω_{x}} f_{r} (x, ω^{'}, ω) \frac{d^{2} Φ_{i} (x, ω^{'})}{d ω_{i}^{'} d A_{i}} d ω_{i}^{'}$

$L_{r} (x, ω) = \int_{Ω_{x}} f_{r} (x, ω^{'}, ω) \frac{d^{2} Φ_{i} (x, ω^{'})}{d A_{i}}$

也可以写成

$L_{r} (x, ω) = \int_{Ω_{x}} f_{r} (x, ω^{'}, ω) d E_{i} (x, ω^{'})$

下面做近似。
输入功率Φi用距离点x最近的n个光子来计算。每个光子的功率是ΔΦp(ωp)。
得到

$L_{r} (x, ω) \approx \sum_{p = 1}^{n} f_{r} (x, ω_{p}^{'}, ω) \frac{Δ Φ_{p} (x, ω_{p})}{Δ A}$

从点x生成一个球，不断变大直到包含n个光子。用此时的半径算A。如图8所示。

$L_{r} (x, ω) \approx \frac{1}{π r^{2}} \sum_{p = 1}^{n} f_{r} (x, ω_{p}^{'}, ω) Δ Φ_{p} (x, ω_{p})$

这个近似包含很多假设，准确度依赖于photon map和这个公式中采用的光子数量。

数量越多越准确。

球体在物体边缘或墙角可能包含错误的光子，造成结果错误。
这些错误区域也是跟光子数量相关。
如果n趋近于无穷大，那么近似可以写成等号。

可以用其他形状代替球体，各有好坏。球体比较好算。

从我测试来看，例如发射20000个光子，生成100000个左右的diffuse光子，4000左右的caustic光子。
4000其实挺少的。计算能量的时候某个交点周围会取n个caustic光子，这个n如何取值？
如果定的比较大，那么会到较大范围去找齐光子，那么r可能会很大，就会有严重问题。
跑老远找光子，贡献又不大，但会把r拖大，造成把整体能量变的很小。
如果定的很小，那么图像就类似青一块紫一块。
还有边界问题。比如在一个突出的墙角，可能会穿过墙壁找光子。

3.1 Filtering#

如果photon map中的光子数量太少，边缘的radiance就会模糊。

这种artifact有时有益，有时有害。

要减少这种模糊，估算的radiance需要做过滤。给越靠近点x的光子增加越多权重。

收集光子用的是球体，但过滤并不是3d的。收集的光子都看作是在同一微小表面。

cone filter#

$ω_{p c} = 1 - \frac{d_{p}}{k r}$

加上过滤和normal

$L_{r} (x, ω) \approx \frac{\sum_{p = 1}^{n} f_{r} (x, ω_{p}^{'}, ω) Δ Φ_{p} (x, ω_{p}) ω_{p c}}{(1 - 2 / 3 k) π r^{2}}$

$k \geq 1$

Gaussian filter#

$ω_{p g} = α (1 - \frac{1 - e^{- β \frac{d_{p}^{2}}{2 r^{2}}}}{1 - e^{- β}})$

$α = 0.918$

$β = 1.953$

$L_{r} (x, ω) \approx \frac{1}{π r^{2}} \sum_{p = 1}^{n} f_{r} (x, ω_{p}^{'}, ω) Δ Φ_{p} (x, ω_{p}) ω_{p g}$

3.2 rendering#

photon map可以独立使用。有几种方法进行可视化。

暂不考media。

总的输出radiance为自发的radiance加反射出的radiance。

$L_{o} (x, ω) = L_{e} (x, ω) + L_{r} (x, ω)$

反射radiance为半球上所有方向的入射radiance，经过cos项和brdf处理，积分。

$L_{r} (x, ω) = \int_{Ω_{x}} f_{r} (x, ω^{'}, ω) L_{i} (x, ω^{'}) | c o s θ_{i} | d ω_{i}^{'}$

L_{r}

可用蒙特卡洛光线追踪算出。即我们之前实现的路径追踪。

这些方法比较耗时。而光子映射更高效。

brdf分成specular和diffuse两部分

$f_{r} (x, ω^{'}, ω) = f_{r, s} (x, ω^{'}, ω) + f_{r, d} (x, ω^{'}, ω)$

$L_{i}$ 分成三部分

$L_{i, l} (x, ω^{'})$ 光源的直接光

$L_{i, c} (x, ω^{'})$ caustics。间接光。至少经过一次specular反射或折射，再打到diffuse面。

$L_{i, d} (x, ω^{'})$ 间接光。至少一次diffuse反射。

$L_{i} (x, ω^{'}) = L_{i, l} (x, ω^{'}) + L_{i, c} (x, ω^{'}) + L_{i, d} (x, ω^{'})$

最终的 $L_{r}$ 由4项组成

$L_{r} (x, ω) = \int_{Ω_{x}} f_{r} (x, ω^{'}, ω) L_{i} (x, ω^{'}) | c o s θ_{i} | d ω_{i}^{'} =$

$\int_{Ω_{x}} f_{r} (x, ω^{'}, ω) L_{i, l} (x, ω^{'}) | c o s θ_{i} | d ω_{i}^{'} +$

$\int_{Ω_{x}} f_{r, s} (x, ω^{'}, ω) (L_{i, c} (x, ω^{'}) + L_{i, d} (x, ω^{'})) | c o s θ_{i} | d ω_{i}^{'} +$

$\int_{Ω_{x}} f_{r, d} (x, ω^{'}, ω) L_{i, c} (x, ω^{'}) | c o s θ_{i} | d ω_{i}^{'} +$

$\int_{Ω_{x}} f_{r, d} (x, ω^{'}, ω) L_{i, d} (x, ω^{'}) | c o s θ_{i} | d ω_{i}^{'}$

其中

直接光#

$\int_{Ω_{x}} f_{r} (x, ω^{'}, ω) L_{i, l} (x, ω^{'}) | c o s θ_{i} | d ω_{i}^{'}$

最重要的部分。

准确计算

用蒙特卡洛光线追踪。如果交点和光源之间有其他物体，那么直接剔除。

另有一种shadow photons方法。
近似计算

用global photon map计算。

Specular和glossy反射#

$\int_{Ω_{x}} f_{r, s} (x, ω^{'}, ω) (L_{i, c} (x, ω^{'}) + L_{i, d} (x, ω^{'})) | c o s θ_{i} | d ω_{i}^{'}$

这里不必使用任何photon map。因为这个值由 $f_{r, s}$ 主导，只在镜面反射方向的狭小范围内起作用。

也可以用pm优化，但需要大量光子。不合算。

Caustics#

$\int_{Ω_{x}} f_{r, d} (x, ω^{'}, ω) L_{i, c} (x, ω^{'}) | c o s θ_{i} | d ω_{i}^{'}$

准确计算

需要从caustics photon map计算功率。光子数量越大结果越准。
近似计算

从global photon map计算。

Multiple Diffuse Reflections#

$\int_{Ω_{x}} f_{r, d} (x, ω^{'}, ω) L_{i, d} (x, ω^{'}) | c o s θ_{i} | d ω_{i}^{'}$

计算至少有一次diffuse反射的光路。因为有diffuse，所以这个值比较“柔和”。

准确计算

要用蒙特卡洛光线追踪。

通常蒙特卡洛计算diffuse间接光比较低效，但是这个积分有一些特殊性质，最终效果还行。

间接光很柔和，因为我们已经把caustics分出去了。caustics通常是蒙特卡洛光追中noise的主因。
近似计算

用global photon map。

可以看到所有近似计算都可以用global photon map。global photon map本质上包含了所有光照！

4. Participating Media#

之前的讨论中，都是考虑光子与mesh表面碰撞，严格来说这只发生在真空中。
即使是很纯净的空气也会改变传播方向。
在更复杂或要求更高的场景中，就不能忽略这个效果，比如烟雾、云等等。
这类介质叫做participating media。

另一类participating media是有些透明度的材质，比如大理石、皮肤、植物等。

涉及到subsurface。要解决这些问题，得解决基本的participating media。

光子映射很适合解决这类问题。是第一个完全模拟出subsurface的算法。

4.1 Light Scattering in Participating Media#

当一个光子进入一种介质，可以不受影响继续往前走，或者在某个位置与介质交互。

交互的结果：要么被吸收，要么被发散出去(scattered)。
发散因子σs(scattering coefficient)
吸收因子σa(absorption coefficient)

这些因子的单位是百分比每单位距离。比如选定0.2/米，可以理解为光子每移动1米，有20%概率被吸收。
或者说100个光子移动1米，有20个会被吸收。
或者说光子在这个介质中，被吸收前，平均能移动1/0.2=5米。

当光穿过介质，可以看作是radiance的连续改变。

在某一点x处，

ω

方向

对于发散(out-scattering)，改变为

$(ω \cdot \nabla) L (x, ω) = - σ_{s} (x) L (x, ω)$

对于吸收，改变为

$(ω \cdot \nabla) L (x, ω) = - σ_{a} (x) L (x, ω)$

合起来，对于 $σ_{t} = σ_{s} + σ_{a}$ 有

$(ω \cdot \nabla) L (x, ω) = - σ_{t} (x) L (x, ω)$

穿过介质的过程中，存在光的in-scattering，所以会有个gain。

$(ω \cdot \nabla) L (x, ω) = σ_{s} (x) \int_{Ω_{4} π} p (x, ω^{'}, ω) L i (x, ω^{'}) d ω^{'}$

ω^{'}

是入射方向。跟渲染方程类似，积分是求所有方向的入射光的贡献。

其中p是phase function，描述光的分布。

另外如果介质本身会发光，又有一个gain。

$(ω \cdot \nabla) L (x, ω) = σ_{a} (x) L_{e} (x, ω)$

一共4个改变，两个加两个减。合起来有

$(ω \cdot \nabla) L (x, ω) =$

$d L (x, ω) = σ_{a} (x) L_{e} (x, ω) + σ_{s} (x) \int_{Ω_{4} π} p (x, ω^{'}, ω) L i (x, ω^{'}) d ω^{'} - σ_{t} (x) L (x, ω)$

4.2 体渲染方程#

要解出

L (x, ω)

有一堆数学上的推导，对于我们数学不熟的人比较跳跃。

论文里直接就给答案，没有推导，只说了一句，对两边同时积分一段距离s。

首先对于 dL(x,ω)=−σa(x)L(x,ω)
L在某个x处的改变量，等于该处的L值乘一个系数。
本质是一个微分方程。我们先要学一下基本的微分方程，一般微积分书里都有。
比如托马斯微积分第16章。
另我的文档一阶微分方程基础学习。其中也包含了一些微积分的重要但容易忽略的知识。

L的这个问题本质和微积分教材里的人口问题、阻力问题是一样的。

都是变化率与实时值相关，然后包含初始值问题。

对于人口问题的微分方程

$\frac{d P}{d t} = k P ， P (0) = P_{0}$

$P = P_{0} e^{k t}$

变量是时长t。但我们的L变量是任意点x，而且形式不一样，能否转成一个标准的微分方程？

对于人口来说是t从0开始，经过一段时间之后P的情况。
类似可以用距离s的思想来重写，表示成光子从0开始
移动一段距离s之后L的情况。
那么

$d L (x, ω) = - σ_{a} (x) L (x, ω)$

转变成微分方程

$\frac{d L (s)}{d s} = - σ_{a} (s) L (s)$

这样就舒服了，可以用这个思想重写整个式子。

$\frac{d L (s)}{d s} = σ_{a} (s) L_{e} (s) + σ_{s} (s) \int_{Ω_{4} π} p (s, ω^{'}) L i (s, ω^{'}) d ω^{'} - σ_{t} (s) L (s)$

写成微分方程标准形式

$L^{'} (s) + σ_{t} (s) L (s) = σ_{a} (s) L_{e} (s) + σ_{s} (s) \int_{Ω_{4} π} p (s, ω^{'}) L i (s, ω^{'}) d ω^{'}$

即

$L^{'} (s) + p (s) L (s) = q (s)$

此时p和q可以说都是确定的，只需解微分方程求L。可参考上述的一阶微分方程基础学习

只包含吸收的情况#

如果简化一下，只考虑吸收，那么q为0。有

$L^{'} (s) + p (s) L (s) = 0$

$L (s) = C e^{- P (s)}$

给定初始值L(0)，可以算出C，进而算出L(s)。

$L (0) = C e^{- P (0)}$

$C = L (0) e^{P (0)}$

$L (s) = L (0) e^{P (0)} e^{- P (s)} = L (0) e^{P (0) - P (s)}$

也可以写成

$L (s) = L (0) e^{\int_{s}^{0} p (t) d t}$ (t只是一个变量符号，可以随意换，写成t避免冲突)

甚至可以不写t

$L (s) = L (0) e^{\int_{s}^{0} p}$

写成从0积到s更直观，加个负号。

$L (s) = L (0) e^{- \int_{0}^{s} p_{s} (t) d t}$

这里的s是经过的距离，即从初始位置积到某个距离s。

如果取L(0)为1，L(s)就是积到s时，对于初始值的倍数！

如果是均匀的介质，Ps为常数k，有

$L (s) = e^{- k s}$

也可以按具体的点考虑，那么就是从点x1积到点x2。

$L (s) = L (x_{1}) e^{- \int_{x_{1}}^{x_{2}} p_{x} (t) d t}$

也就是上一节说的两种思想，一个是针对距离，一个是针对某个坐标。

那么p相应地应该有两个版本

p_{s}

和

p_{x}

。最后的结果肯定一样。

完整的情况#

$L^{'} (s) + p (s) L (s) = q (s)$

$L (s) = \int_{a}^{s} e^{- \int_{t}^{s} p (t^{'}) d t^{'}} q (t) d t + C e^{- \int_{a}^{s} p (t^{'}) d t^{'}}$

按[一阶微分方程学习]最后的讨论。选定a=0，结合初始值L(s=0)确定C

$L (0) = C$

代回原方程得到最终的体渲染方程

$L (s) = \int_{0}^{s} e^{- \int_{t}^{s} p (t^{'}) d t^{'}} q (t) d t + L (0) e^{- \int_{0}^{s} p (t^{'}) d t^{'}}$

意思就是一个初始能量为L(0)的光子，在介质中移动s距离后的最终能量。

4.3 Phase Function#

phase function描述散射的光如何在介质中分布。

p定义再某点上某个方向w对某个方向w’的贡献值。

p必须满足给定一个方向w，所有方向w’射过来的能量对w方向的贡献，积起来值为1。

$\int_{Ω_{4} π}^{} p (x, ω^{'}, ω) d ω^{'} = 1$

p和brdf相似，但是有两个大不同： p没有单位，且是normalized。

可以用入射光和出射光的夹角 $θ$ 简写

各向同性#

是个常数，各个方向均匀分布

$p (θ) = \frac{1}{4 π}$

Henyey-Greenstein#

最常用的phase function。最初被用来解释星尘的散射，后来被用来描述海洋、云、皮肤、石头等等。

$p (θ) = \frac{1 - g^{2}}{4 π (1 + g^{2} - 2 g c o s θ)^{1.5}}$

正的g为前向散射，负的g为后向散射，0为各向同性。

θ = 0

为前向，

θ = π

为反向。

可以画出不同g值对应的图形。可以看到g=0时，是一个圆形，任何角度都是1/4π。
g2越趋近于1，圆形越扁，越分布在正向和反向附近。
g值正地越多越往正向集中，负地越多越往负向集中。

4.4 Ray Marching#

体渲染方程，除了一些简单场景之外，只能用数值积分来计算。
L是关于距离s的函数，很自然想到模拟光子不断前进，就是所谓的Ray Marching。
把前进路线分成很多小步，每步算一下体渲染方程。必要时做一些简化。

那么就是一个近似的解，细分越多越精确。

啊照这样说，那是不是根本就不用推导那一大片微分方程？每一步直接算最初的dL和更新L不就行了？

4.4.1 自适应ray matching#

上面的ray marching时平均步长。如果是非均匀的材质，或者有本地波动(比如阴影)，最好用自适应步长。例如当材质特性发生较大变化时，可采用较小的步长，以更好地反映出细节。

方法1 可以以递归的形式来操作，从两个端点开始，从中间破开，左右递归下去。

直到距离达到设定的限制，或者两个端点的差别小到一定程度(不用继续破开)。

方法2 根据sigma_t来调整步长。这个对于未知的材质可能很有效。

取均匀分布的随机数

貌似涉及均匀采样

4.5 Photon Tracing#

4.5.1 Photon Scattering#

两种情况光子与介质交互，或发散或被吸收。光子与某表面交互，反射或传输或吸收。

光子进入介质时，不看作与表面交互，而是直接往里面走，进行与介质的交互。

之前已提到sigma的意义，光子进行下一个交互前平均移动的距离。

光子如果发生scattering，对phase函数进行采样，得到新的方向。

推导见(从0开始学习PBRT)。

4.5.2 photon storing#

每当光子与介质进行交互，就要存光子，无论是吸收还是散射。

对于介质，我们起一个独立的volume photon map。因为介质交互中的radiance的计算方式与表面交互不一样。

一个优化是，只存至少散射过一次的光子。这样省去了直接光的部分。直接光可以用传统ray tracing高效算出，不需要费劲用光子来算。

4.5.3 photon emission#

光子也可以从介质中发出。

4.6 radiance的计算#

4.7 渲染#

4.8 次表面散射#

次表面散射存在于所有非金属材质。

一般会简化掉，做成光子碰到表面后立即从交点往外反射离开，而不是往里面散射。这样简化对于皮肤、玉石、牛奶等材质来说，效果就很差。

要准确地处理translucent材质，就需要遵循事实：光从某点进入，会经过次表面散射，从其他点出去。

photon mapping是第一种完全实现这个原理的算法。

photon tracing#

同上述的photon tracing

渲染#

结果#

还有较多的细节没仔细看。

目前只做了点光源。

光子能量不知道怎么取值。大概加了个scale，比如5000，能得到还算正常的图。

difffuse光子数量通常是caustic数量的5倍。

光子数量问题，较难调整。数量少的话随机性太大。多的话很耗时。

等换了新电脑再好好试试。

knn中k的取值，是取周围5个光子，10个还是100个500个？
没有明确的说法。一般是说越大越准。但有一系列问题。
比如墙角，如果光子密度小，又要强行取一批光子，那么可能越过墙角，明显是错的。
也就是光子取的太多，范围可能太大，中间可能有遮挡。
取的少么，又可能不够准。
最好就是发射巨量的光子，使得每个交点周围的很小范围内都能找到需要的光子。

它这个渲染本质上是看光子的密度。那么限定k和限定r其实是一回事？

资料#

[Realistic Image Synthesis Using Photon Mapping]

https://www.scratchapixel.com/lessons/3d-basic-rendering/volume-rendering-for-developers/volume-rendering-summary-equations

https://www.cs.cornell.edu/courses/cs667/2005sp/notes/08tchao.pdf
https://www.cs.cornell.edu/courses/cs667/2005sp/notes/09wang.pdf
这个资料总的很有帮助，有较多推导。但是有几个符号和字母小错，而且有的推导貌似不严密。

论文学习#

`Progressive Photon Mapping <http://graphics.ucsd.edu/~henrik/papers/progressive_photon_mapping/progressive_photon_mapping.pdf>`__

基本光子映射的一个问题是需要大量的光子来达到较高的精确度，相应地需要大量的memory。

ppm是一个多pass算法，可以不断迭代以增加准确度，同时不必维护每次迭代的photon map。从而不必无限制地增加memory用量。

基本光子映射的 $L_{r}$

$L_{r} (x, ω) \approx \frac{1}{π r^{2}} \sum_{p = 1}^{n} f_{r} (x, ω_{p}^{'}, ω) Δ Φ_{p} (x, ω_{p})$

这只是一个近似。只有当n趋近于无穷大时，才会收敛到一个正确的无偏结果。

ppm算法第一遍是ray tracing，后续每个pass都是photon tracing。

每次photo tracing都增加了gi的准确度。所以称为progressive。

ray tracing pass#

和蒙特卡洛光追类似，从眼睛向每个pixel发光，直到碰到第一个非specular面。

对于每个hitpoint，存下列信息。

交点位置
normal
入射光方向
brdf
pixel坐标
pixel权重
当前统计半径
光子累积数量
累积能量

photon tracing pass#

任务是对ray tracing的结果数据中的能量进行累积。

每一遍可以追踪一定数量的光子，建立photon map。

每一遍结束时，遍历所有交点，找到交点周围的光子。
用这些光子更新能量。
这样形成一个不断完善gi的过程。直到图像达到满意的质量。

上一遍用过的光子以后就不需要了，可以扔掉。

每一遍过后都可以更新图像。

能量的估算#

基本光子映射的 $L_{r}$

$L_{r} (x, ω) \approx \frac{1}{π r^{2}} \sum_{p = 1}^{n} f_{r} (x, ω_{p}^{'}, ω) Δ Φ_{p} (x, ω_{p})$

依赖于光子的密度

$d (x) = \frac{n}{π r^{2}}$

其中假定n个光子分布在交点所处的微小disk面上，半径为r。

如果再生成一个新的map，又可以算一个新的密度

$d_{n e w} (x) = \frac{n_{n e w}}{π r^{2}}$

r不变。两个密度平均一下能得到更准确的密度。
要注意这个方法并不是consistent的，且不会准确收敛在x上。
它只是在r上算平均值。没法体现r之内的细节。且准确度依赖于每一遍的map。

按这个方法处理一系列photon map就会收敛到正确结果。

其中的一个关键点是，减小估算半径，同时累加光子。

半径消减#

给x定义一个半径 $R (x)$ 。目标是增加其中的光子数 $N (x)$ 同时减小半径。

假设已有 $N (x)$ 个光子。此时经过一个pass又找到 $M (x)$ 个光子。可计算新密度

$d^{'} (x) = \frac{N (x) + M (x)}{π R (x)^{2}}$

接着要缩小半径

$R^{'} (x) = R (x) - R_{c u t}$

假设 $R (x)$ 中的光子密度处处相等。那么新的半径中的光子数量为

$N^{'} (x) = π R^{'} (x)^{2} d^{'} (x) = π (R (x) - R_{c u t})^{2} d^{'} (x)$

另一方面，为了结果能收敛，每一遍我们必须增加光子。
引入α来控制实际的加入量。可以随便自定这个量，并不是定死M。
M是计算密度改变这个事实，而α是用于估算。

$N^{'} (x) = N (x) + α M (x)$

几个方程联立可解出 $R_{c u t}$ 和 $R^{'} (x)$

$π (R (x) - R_{c u t})^{2} \frac{N (x) + M (x)}{π R (x)^{2}} = N (x) + α M (x)$

$R (x) - R_{c u t} = R (x) \sqrt{\frac{N (x) + α M (x)}{N (x) + M (x)}} = R^{'} (x)$

这样就可以用

α

来控制每次的光子添加量和半径的缩减量。

直接指定半径缩减量，反推

α

也是可以的。用

α

更方便。

注意新的半径 $R^{'} (x)$ 对于每个交点都要算。

功率的矫正#

设某交点现有光子的功率为

F_{N}

，一遍过后新增功率

F_{M}

。

做半径缩减后，新的功率简单按面积的比例计算即可

$F_{N}^{'} = F l u x_{n e w} = (F_{N} + F_{M}) \frac{π R^{'} (x)^{2}}{π R (x)^{2}}$

$= (F_{N} + F_{M}) \frac{(R (x) \sqrt{\frac{N (x) + α M (x)}{N (x) + M (x)}})^{2}}{R (x)^{2}}$

$= (F_{N} + F_{M}) \frac{N (x) + α M (x)}{N (x) + M (x)}$

辐射率估算#

基本光子映射方程

$L_{r} (x, ω) \approx \frac{1}{π r^{2}} \sum_{p = 1}^{n} f_{r} (x, ω_{p}^{'}, ω) Δ Φ_{p} (x, ω_{p})$

其中求和即与 $F$ 等价。

$L_{r} (x, ω) \approx \frac{1}{π R (x)^{2}} F_{N} (x)$

可以看出本质就是把执行一次的大量光子计算(pm)，转化为了一遍遍增加光子，同时增强画质(ppm)。

最后要normalize

$L_{r} (x, ω) \approx \frac{F_{N} (x)}{π R (x)^{2} \times 总光子数}$

总结#

与pm有不少区别。

每个hitpoint设定一个初始r值。
每次添加光子，比如2000个。计算每个hitpoint在r范围内新增的光子数。
如果这样那么knn就是不必要的？因为必须算出r范围内所有新增光子，而不是自定义取最近的k个光子。

论文学习#

`Stochastic Progressive Photon Mapping <http://graphics.ucsd.edu/~henrik/papers/sppm/stochastic_progressive_photon_mapping.pdf>`__

ppm#

ppm是一个多pass算法，用统计的方法累积光子，来解决渲染方程。

radiance可以写成

$L_{r} (x, ω) \approx \frac{F_{i} (x)}{N_{e} (i) π R_{i} (x)^{2}}$

其中的i为当前pass数，每个pass中F，N，R都得到更新，结果变得更好。

当i趋近于无穷大时L达到准确值(收敛)。

$L (x, ω) = lim_{i \to \infty} \frac{F_{i} (x)}{N_{e} (i) π R_{i} (x)^{2}}$

sppm#

sppm的动机是，计算一块区域内的准确平均radiance。

比如motion blur需要算平均值。depth-of-field也需要算平均值。