目录
Basic Radiometry
整理一下渲染方程相关的基础知识.
Radiant Energy and Flux (Power)
-
定义:Radiant Energy (电磁辐射的能量)。以焦耳为单位,符号表示: Q [ J = J o u l e ] Q[J=Joule] Q[J=Joule]
-
定义:Flux(单位时间的能量): 单位时间内发射、反射、传输或接收的能量。
Φ ≡ d Q d t [ W = W a t t ( 功 率 ) ] [ l m = l u m e n ( 流 明 ) ] ∗ (1) \tag1\Phi\equiv\frac{dQ}{dt}[W=Watt(功率)][lm=lumen(流明)]^* Φ≡dtdQ[W=Watt(功率)][lm=lumen(流明)]∗(1)- 光学中使用流明 l m lm lm代表功率 W W W
- ≡ \equiv ≡ 恒等
- Flux还可以解释为单位时间内通过某感光介质光子的数量
-
辐射能量 方向性的能量单位 Radiant Intensity or Intensity
- power per unit solid angle
- 单位立体角的能量
- I ( ω ) ≡ d Φ d ω I(\omega)\equiv\frac{d\Phi}{d\omega} I(ω)≡dωdΦ
-
I
(
ω
)
≡
d
Φ
d
ω
[
W
s
r
]
[
l
m
s
r
=
c
d
=
c
a
n
d
e
l
a
]
I(\omega)\equiv\frac{d\Phi}{d\omega}\begin{bmatrix}\frac{W}{sr}\end{bmatrix}\begin{bmatrix}\frac{lm}{sr}=cd=candela\end{bmatrix}
I(ω)≡dωdΦ[srW][srlm=cd=candela]
- Φ \Phi Φ是能量
- ω \omega ω是立体角
- l m lm lm为流明
- candela是七种国际单位制基本单位之一
-
立体角
- 球面的面积除以半径的平方
Ω
=
A
r
2
\Omega = \frac{A}{r^2}
Ω=r2A
- A A A球面中某一区域的面积
- r r r为球的半径
- 整个球的立体角 4 π 4\pi 4π
- 立体角单位steradians
- 球面的面积除以半径的平方
Ω
=
A
r
2
\Omega = \frac{A}{r^2}
Ω=r2A
-
立体角的定义
- 单位立体角或微分立体角
- 单位球面面积
-
d
A
=
(
r
d
θ
)
(
r
s
i
n
θ
d
ϕ
)
=
r
2
s
i
n
θ
d
θ
d
ϕ
dA = (rd\theta)(rsin\theta d\phi) = r^2sin\theta d\theta d\phi
dA=(rdθ)(rsinθdϕ)=r2sinθdθdϕ
-
θ
\theta
θ是方向向量与Z轴的夹角
-
θ
\theta
θ的求解
- 将球上一点投影到ZY平面上, 投影点与原点形成的向量与Z轴的夹角
-
θ
\theta
θ的求解
- r d θ rd\theta rdθ: 以球心为圆心, r r r为半径, Z轴所在的平面(任意平面)中的弧长
-
ϕ
\phi
ϕ
- 将球上的一点投影到XY平面上, 投影点与原点形成的向量与X轴的夹角
- r s i n θ d ϕ rsin\theta d\phi rsinθdϕ: 球上任意一点投影到Z轴上, 以该投影点为圆心,以 r s i n θ rsin\theta rsinθ为半径, 垂直于Z轴的平面上圆的弧长
- d A dA dA: 球面上的单位面积, 就是两段弧长的乘积
-
θ
\theta
θ是方向向量与Z轴的夹角
-
d
A
=
(
r
d
θ
)
(
r
s
i
n
θ
d
ϕ
)
=
r
2
s
i
n
θ
d
θ
d
ϕ
dA = (rd\theta)(rsin\theta d\phi) = r^2sin\theta d\theta d\phi
dA=(rdθ)(rsinθdϕ)=r2sinθdθdϕ
- 单位立体角
- d ω = d A r 2 = s i n θ d θ d ϕ d\omega=\frac{dA}{r^2} = sin\theta d\theta d\phi dω=r2dA=sinθdθdϕ
- 整个球的立体角
- Ω = ∫ S 2 d ω = ∫ 0 2 π ∫ 0 π s i n θ d θ d ϕ = 4 π \Omega = \int_{S^2}d\omega = \int_0^{2\pi}\int_0^\pi sin\theta d\theta d\phi = 4\pi Ω=∫S2dω=∫02π∫0πsinθdθdϕ=4π
-
辐射度量学中使用 ω \omega ω表示三维空间中的方向(单位向量)
-
Intensity
- 光源在任何一个方向上的亮度
- 所有方向上的亮度除以对应单位立体角
- 所有单位立体角的Intensity积分起来即可得到Flux(
Φ
\Phi
Φ)
- Φ = ∫ S 2 I d ω = 4 π I \Phi=\int_{S^2}Id\omega = 4\pi I Φ=∫S2Idω=4πI
- 设一个点光源均匀的辐射能量, 任何一个方向上的Intensity就是 I = Φ 4 π I = \frac{\Phi}{4\pi} I=4πΦ
- 光源在任何一个方向上的亮度
-
设灯泡流明为815, 求Intensity
- I n t e n s i t y = 815 ( l u m e n s ) / 4 π ( s r ) = 65 c a n d e l a s Intensity = 815(lumens)/4\pi(sr) = 65 candelas Intensity=815(lumens)/4π(sr)=65candelas
-
Irradiance: power per unit area
- 物体表面接受多少光的能量 (各个方向辐射)
-
E
(
x
)
=
d
Φ
(
x
)
d
A
E(x) = \frac{d\Phi(x)}{dA}
E(x)=dAdΦ(x)
- Φ ( x ) \Phi(x) Φ(x)为能量, A A A为物体表面的面积
- irradiance表示一个极小的能量除以一个极小的面积所得到的能量.
- [ W m 2 ] [ l m m 2 = l u x ] \begin{bmatrix}\frac{W}{m^2}\end{bmatrix}\begin{bmatrix}\frac{lm}{m^2}=lux\end{bmatrix} [m2W][m2lm=lux]
-
Irradiance与Intensity区别
- Intensity是单位立体角包含的能量
- Irradiance是单位面积所包含的能量
- 该面积是指光源投影到垂直方向上的面积, 因此光源垂直于平面能量获取的最多
- E = Φ A c o s θ E = \frac{\Phi}{A}cos\theta E=AΦcosθ
- Intensity在每个传播方向上的能量是不会发生变化的
- Irradiance在传播过程会衰减的
- 当光传播单位距离时(1), irradiance为:
- E = Φ 4 π ∗ 1 E = \frac{\Phi}{4\pi*1} E=4π∗1Φ
- 当光线传播距离为
r
r
r时, irradiance为:
- E = Φ 4 π ∗ r 2 E = \frac{\Phi}{4\pi*r^2} E=4π∗r2Φ
- 可以看出,光线在传播过程中,其所对于的球面面积逐渐增大,因此irradiance随着光线传播会越来越小。
- 当光传播单位距离时(1), irradiance为:
-
Radiance:power per unit solid angle, per porjectged unit area
- 是光在传播中能量的度量
- 描述光线属性
L ( p , ω ) ≡ d 2 Φ ( p , ω ) d ω d A c o s θ (2) \tag2L(p, \omega) \equiv \frac{d^2\Phi(p,\omega)}{d\omega dAcos\theta} L(p,ω)≡dωdAcosθd2Φ(p,ω)(2)- d A c o s θ dAcos\theta dAcosθ 投影面积
- [ W s r m 2 ] = [ c d m 2 = l m s r m 2 = n i t ] \begin{bmatrix}\frac{W}{srm^2}\end{bmatrix} = \begin{bmatrix}\frac{cd}{m^2} = \frac{lm}{srm^2} = nit\end{bmatrix} [srm2W]=[m2cd=srm2lm=nit]
- 从公式(2)可以看出:
- Radiance可以看作:Irradiance per solid angle
- Radiance可以看作:intensity per projected unit area
- 个人理解, Radiance其实就是光在某一方向上辐射的能量(单位立体角, 单位投影面积, 即 Δ \Delta Δ)
-
Incident Radiance
- 入射光线的能量(单位面积接受)
- L ( p , ω ) = d E ( p ) d ω c o s θ L(p, \omega) = \frac{dE(p)}{d\omega cos\theta} L(p,ω)=dωcosθdE(p)
- 从上述公式可以看出, Incident Radiance 是irradiance 在单位立体角中获取的能力
- c o s θ cos\theta cosθ表示,入射光线与单位面积平面的夹角。因为 d A dA dA接受的能量与入射角度相关
- 入射光线的能量(单位面积接受)
-
Exiting Randiance
-
d
A
dA
dA往某个方向发射出的能力
- 离开单位投影面积平面的intensity
- L ( p , ω ) = d I ( p , ω ) d A c o s θ L(p, \omega) = \frac{dI(p, \omega)}{dAcos\theta} L(p,ω)=dAcosθdI(p,ω)
- 从上述公式可以看出, Exiting Radiance 是irradiance 在某个方向上的能量
- c o s θ cos\theta cosθ表示,出射光线与单位面积平面的夹角。
-
d
A
dA
dA往某个方向发射出的能力
-
Irradiance 与 Radiance
- Irradiance: d A dA dA收到的所有能量
- Radiance: 某一方向接收到的能量
- 因此, Irradiance可以看做是各个方向的Radiance积分起来能量
d E ( p , ω ) = L i ( p , ω ) c o s θ d ω (3-1) \tag{3-1}dE(p,\omega) = L_i(p, \omega) cos\theta d\omega dE(p,ω)=Li(p,ω)cosθdω(3-1)
E ( p ) = ∫ H 2 L i ( p , ω ) c o s θ d ω (3-2) \tag{3-2}E(p) = \int_{H^2}L_i(p,\omega)cos\theta d\omega E(p)=∫H2Li(p,ω)cosθdω(3-2)- L i ( p , ω ) L_i(p, \omega) Li(p,ω) 为在 ω \omega ω方向上的Radiance
- H 2 H^2 H2表示半球(光传播只考虑上半球辐射的能量)
- 公式(3)的意义是解释Irradiance和Radiance之间的关系
Bidirectional Reflectance Distribution Function (BRDF)
- BRDF: 双向反射分布函数
- 意义
- 光从某一方向进入物体表面, 有多少能量进入
- 反射到某个方向上的能量是多少
- 不同的反射方向会分布多少能量
- 意义
- Reflection at a point
- 解释: 从 ω i \omega_i ωi来的Radiance 打到 d A dA dA上, d A dA dA吸收能量, 并转化发射出能量到另一个方向上
- 公式解释
-
d
E
(
ω
i
)
=
L
(
ω
i
)
c
o
s
θ
i
d
ω
i
dE(\omega_i) = L(\omega_i)cos\theta_id\omega_i
dE(ωi)=L(ωi)cosθidωi
- ω i \omega_i ωi方向上接受的能量等于所有 L ( ω i ) L(\omega_i) L(ωi)(Radiance)总和
-
d
L
r
(
ω
r
)
dL_r(\omega_r)
dLr(ωr)
- ω i \omega_i ωi方向接受所有的能量 d E ( ω i ) dE(\omega_i) dE(ωi)反射出去的能量为: 在 ω r \omega_r ωr方向的Irradiance
-
d
E
(
ω
i
)
=
L
(
ω
i
)
c
o
s
θ
i
d
ω
i
dE(\omega_i) = L(\omega_i)cos\theta_id\omega_i
dE(ωi)=L(ωi)cosθidωi
- BRDF
- 是
d
A
dA
dA收到的能量(Irradiance), 反射到另一个方向上的能量
f r ( ω i → ω r ) = d L r ( ω r ) d E i ( ω i ) = d L r ( ω r ) L i ( ω i ) c o s θ i d ω i (4) \tag4 f_r(\omega_i\to \omega_r) = \frac{dL_r(\omega_r)}{dE_i(\omega_i)} = \frac{dL_r(\omega_r)}{L_i(\omega_i)cos\theta_id\omega_i} fr(ωi→ωr)=dEi(ωi)dLr(ωr)=Li(ωi)cosθidωidLr(ωr)(4)- f r ( ω i → ω r ) f_r(\omega_i\to \omega_r) fr(ωi→ωr)表示单位面积平面(着色点)从 ω i \omega_i ωi接受到的Irradiance反射到 ω r \omega_r ωr方向上的能量占Irradiance的比例
- d L r ( ω r ) dL_r(\omega_r) dLr(ωr)为反射到 ω r \omega_r ωr方向上的能量. 即 ω r \omega_r ωr方向的Radiance的微分
- L i ( ω i ) c o s θ i d ω i L_i(\omega_i)cos\theta_id\omega_i Li(ωi)cosθidωi为 ω i \omega_i ωi方向上接受所有Radiance的总和( i i i为上半球任意方向)
- 是
d
A
dA
dA收到的能量(Irradiance), 反射到另一个方向上的能量
The Reflection Equation
- 反射方程
- 任意着色点在各种光照下, 任何一个光线进入方向对出射光线的贡献, 所有入射光贡献加起来, 既是最终反射出的能量
- 反射方程不仅反射各个方向上光源的能量(Radiance), 也可反射其他面片反射出的能量
L r ( p , ω r ) = ∫ H 2 f r ( p , ω i → ω r ) L i ( p , ω i ) c o s θ i d ω i (5) \tag5 L_r(p, \omega_r) = \int_{H^2}f_r(p, \omega_i \to\omega_r)L_i(p, \omega_i)cos\theta_id\omega_i Lr(p,ωr)=∫H2fr(p,ωi→ωr)Li(p,ωi)cosθidωi(5)- L r ( p , ω r ) L_r(p, \omega_r) Lr(p,ωr)为 ω r \omega_r ωr方向上反射出的Radiance. 也就是着色点在所有可能的入射光下, 最后反射出去的能量
- H 2 H^2 H2表示上半球所有方向(角度)
- f r ( p , ω i → ω r ) f_r(p, \omega_i \to\omega_r) fr(p,ωi→ωr)表示从 ω i \omega_i ωi入射方向反射到 ω r \omega_r ωr方向上的BRDF
- L i ( p , ω i ) L_i(p, \omega_i) Li(p,ωi)表示着色点在 ω i \omega_i ωi入射方向接收到的Radiance
- d ω i d\omega_i dωi表示上半球上每个方向
The Rendering Equation
L o ( p , ω o ) = L e ( p , ω o ) + ∫ Ω + L i ( p , ω i ) f r ( p , ω i , ω o ) ( n ⋅ ω i ) d ω i (6) \tag6 L_o(p, \omega_o) =L_e(p, \omega_o) + \int_{\Omega^+}L_i(p, \omega_i)f_r(p, \omega_i, \omega_o)(n\cdot\omega_i)d\omega_i Lo(p,ωo)=Le(p,ωo)+∫Ω+Li(p,ωi)fr(p,ωi,ωo)(n⋅ωi)dωi(6)
- 其中
- L o ( p , ω o ) L_o(p, \omega_o) Lo(p,ωo) 表示着色点 p p p在 ω o \omega_o ωo方向上最终反射出的能量
- L e ( p , ω o ) L_e(p, \omega_o) Le(p,ωo)表示自发光辐射出的能量
- Ω + \Omega^+ Ω+与 H 2 H^2 H2相同, 表示上半球的每个方向
- L i ( p , ω i ) L_i(p, \omega_i) Li(p,ωi)表示着色点在 ω i \omega_i ωi入射方向接收到的Radiance
- f r ( p , ω i , ω o ) f_r(p, \omega_i, \omega_o) fr(p,ωi,ωo)表示表示从 ω i \omega_i ωi入射方向反射到 ω o \omega_o ωo方向上的BRDF
- ( n ⋅ ω i ) (n\cdot\omega_i) (n⋅ωi)表示角度, n n n为法线, ω i \omega_i ωi表示入射方向
- 全局光照: 直接光照和间接光照的集合
- 渲染方程的解释
- 单个点光源
- L r ( x , ω r ) = L e ( x , ω r ) + L i ( x , ω i ) f r ( x , ω i , ω r ) ( ω i , n ) L_r(x, \omega_r) =L_e(x, \omega_r) + L_i(x, \omega_i)f_r(x, \omega_i, \omega_r)(\omega_i,n) Lr(x,ωr)=Le(x,ωr)+Li(x,ωi)fr(x,ωi,ωr)(ωi,n)
- 着色点 x x x在 ω r \omega_r ωr方向上最终反射出的能量等于自发光能量加上着色点接收点光源的能量
- 多个点光源
- L r ( x , ω r ) = L e ( x , ω r ) + ∑ L i ( x , ω i ) f r ( x , ω i , ω r ) ( ω i , n ) L_r(x, \omega_r) =L_e(x, \omega_r) + \sum L_i(x, \omega_i)f_r(x, \omega_i, \omega_r)(\omega_i,n) Lr(x,ωr)=Le(x,ωr)+∑Li(x,ωi)fr(x,ωi,ωr)(ωi,n)
- 着色点 x x x在 ω r \omega_r ωr方向上最终反射出的能量等于自发光能量加上着色点接收所有点光源的能量之和
- 面光源
- L r ( x , ω r ) = L e ( x , ω r ) + ∫ Ω L i ( x , ω i ) f r ( x , ω i , ω r ) c o s θ i d ω i L_r(x, \omega_r) =L_e(x, \omega_r) + \int_\Omega L_i(x, \omega_i)f_r(x, \omega_i, \omega_r)cos\theta_id\omega_i Lr(x,ωr)=Le(x,ωr)+∫ΩLi(x,ωi)fr(x,ωi,ωr)cosθidωi
- 着色点 x x x在 ω r \omega_r ωr方向上最终反射出的能量等于自发光能量加上着色点接收面光源的能量的积分(面片所有点看做点光源)
- 其他物体表面反射来的能量
L r ( x , ω r ) = L e ( x , ω r ) + ∫ Ω L r ( x ′ , − ω i ) f r ( x , ω i , ω r ) c o s θ i d ω i (7-1) \tag{7-1}L_r(x, \omega_r) =L_e(x, \omega_r) + \int_\Omega L_r(x', -\omega_i)f_r(x, \omega_i, \omega_r)cos\theta_id\omega_i Lr(x,ωr)=Le(x,ωr)+∫ΩLr(x′,−ωi)fr(x,ωi,ωr)cosθidωi(7-1)-
着色点 x x x在 ω r \omega_r ωr方向上最终反射出的能量等于自发光能量加上从所有其他着色点接收到的能量积分
-
*对公式中 L r ( x ′ , − ω i ) L_r(x', -\omega_i) Lr(x′,−ωi)存在疑问. 为什么是其他物体表面着色点 x ′ x' x′在 − ω i -\omega_i −ωi方向上的Radiance?
- 搜了一些资料, 发现该公式中使用着色点 x ′ x' x′在 − ω i -\omega_i −ωi方向上的Radiance应该是错误的
- 来源: https://wenku.baidu.com/view/bbc49e093b3567ec102d8ae6.html
- 其中Page3 对渲染方程介绍中, 对其他物体表面出射的Radiance有如下公式
- L i ( x , ω ) = L o ( x ∗ ( x , − ω ) , ω ) L_i(x, \omega) = L_o(x^*(x, -\omega), \omega) Li(x,ω)=Lo(x∗(x,−ω),ω)
- 解释
- 对于着色点 x x x在 ω \omega ω入射方向接收其他物体表面反射来的Radiance L i ( x , ω ) L_i(x, \omega) Li(x,ω)有, 当前着色点 x x x在 − ω -\omega −ω方向上Path tracing的着色点 x ∗ x^* x∗. 而最终接收到源物体表面的Radiance为 L o ( x ∗ ( x , − ω ) , ω ) L_o(x^*(x, -\omega), \omega) Lo(x∗(x,−ω),ω).
- 因此公式中 L r ( x ′ , − ω i ) L_r(x', -\omega_i) Lr(x′,−ωi)是错误的, 应该是 L r ( x ′ , ω i ) L_r(x', \omega_i) Lr(x′,ωi), x ′ x' x′为Path tracing的着色点
- 因此最后公式改写为
L r ( x , ω r ) = L e ( x , ω r ) + ∫ Ω L r ( x ′ ( x , − ω i ) , ω i ) f r ( x , ω i , ω r ) c o s θ i d ω i (7-2) \tag{7-2}L_r(x, \omega_r) =L_e(x, \omega_r) + \int_\Omega L_r(x'(x, -\omega_i), \omega_i)f_r(x, \omega_i, \omega_r)cos\theta_id\omega_i Lr(x,ωr)=Le(x,ωr)+∫ΩLr(x′(x,−ωi),ωi)fr(x,ωi,ωr)cosθidωi(7-2)
- 个人理解, 如果这种解释有误, 请忽视该质疑点
- 温故而知新#109
-
公式中 L r ( x ′ , − ω i ) L_r(x', -\omega_i) Lr(x′,−ωi)有另一种解释: 在渲染中通常将入射光线和出射光线的方向定义为从着色点指向进入或反射方向. 只有在这个前提下公式(7-1)是正确的
-
- 单个点光源
- 渲染方程求解:
- 从公式(7)中可以看出
- L r ( x , ω r ) L_r(x, \omega_r) Lr(x,ωr)未知项
- 自发光项 L e ( x , ω r ) L_e(x, \omega_r) Le(x,ωr)已知
- L r ( x ′ ( x , − ω i ) , ω i ) L_r(x'(x, -\omega_i), \omega_i) Lr(x′(x,−ωi),ωi)未知
- f r ( x , ω i , ω r ) f_r(x, \omega_i, \omega_r) fr(x,ωi,ωr)已知(预处理获得)
- c o s θ i cos\theta_i cosθi已知 (入射光线与法线夹角)
- 公式(7)是第二类的Fredholm积分方程(具体解释忽略,不是数学专业,暂不研究)
- 第二类Fredholm(弗雷德霍姆)积分:
- φ ( x ) + λ ∫ b a K ( x , s ) φ ( s ) d s = f ( x ) \varphi(x)+\lambda\int_b^aK(x,s)\varphi(s)ds = f(x) φ(x)+λ∫baK(x,s)φ(s)ds=f(x)
- 其中
- φ ( x ) \varphi(x) φ(x) 是未知函数
- λ \lambda λ是非0参数
- ( a , b ) (a, b) (a,b)为常数区间或无穷区间
- K ( x , s ) K(x,s) K(x,s)是已知连续函数, 称为积分核
- f ( x ) f(x) f(x)是已知函数, 称为自由项
- 渲染方程为第二类Fredholm积分, 因此可以使用行列式方式求解线性方程组
- 将公式(7)简写为
l ( u ) = e ( u ) + ∫ l ( v ) K ( u , v ) d v (8) \tag8 l(u) = e(u) + \int l(v)K(u,v)dv l(u)=e(u)+∫l(v)K(u,v)dv(8) - 其中
- f r ( x , ω i , ω r ) c o s θ i f_r(x, \omega_i, \omega_r)cos\theta_i fr(x,ωi,ωr)cosθi为积分核,(BRDF和夹角已知)
- L e ( x , ω r ) L_e(x, \omega_r) Le(x,ωr)为自由项(着色自发光是已知的)
- 公式(8)可以改写为矩阵方程(老本行)
L = E + K L (9) \tag9L = E + KL L=E+KL(9) - 求解:
- L = E + K L L = E + KL L=E+KL
- ⇒ I L − K L = E \Rarr IL - KL = E ⇒IL−KL=E (I为单位矩阵)
- ⇒ ( I − K ) L = E \Rarr (I - K)L = E ⇒(I−K)L=E
- ⇒ L = ( I − K ) − 1 E \Rarr L = (I - K)^{-1}E ⇒L=(I−K)−1E
- 根据二项式定理(binomial theorem)可得
- ⇒ ( I + K + K 2 + K 3 + . . . ) E \Rarr (I+K+K^2+K^3+...)E ⇒(I+K+K2+K3+...)E
- ⇒ I E + K E + K 2 E + K 3 E + . . . \Rarr IE+KE+K^2E+K^3E+... ⇒IE+KE+K2E+K3E+...
- 因此公式(9)可以近似写成 (即, 渲染方程可以近似写成)
L = E + K E + K 2 E + K 3 E + . . . (10) \tag{10} L= E+KE+K^2E+K^3E+... L=E+KE+K2E+K3E+...(10) - 其中
- E E E为着色点自发光项
- K E KE KE为直接光照
- K 2 E K^2E K2E为间接光照, 即光线经过一次弹射(One bounce)的能量
- K 3 E K^3E K3E为间接光照, 两次弹射(two bounce)的能量
- 从上述推导可以看出. 渲染方程可以使用公式(10)求出近似解
- 从公式(7)中可以看出
Monte Carlo method
- 蒙特卡罗方法的定义有很多, 自己的理解如下:
- 当求某一事件出现的概率或者某一随机变量的期望时, 首先构造描述这一事件的模型, 然后从已知的概率分布中进行抽样. 最后根据所有抽样样本确定问题的解.
- 蒙特卡罗方法的思想就是, 如果要求解某一问题的解时, 可以先建立一个符合该问题的数学模型, 然后从该数学模型中取某一随机变量(问题解). 由于单个解无法诠释该问题, 但按照某种概率分布随机抽样求解的话, 随着样本数量的增多, 所有解(样本)越能诠释问题.从而获得问题的解.
- 数学应用:
- 通常蒙特·卡罗方法通过构造符合一定规则的随机数来解决数学上的各种问题。对于那些由于计算过于复杂而难以得到解析解或者根本没有解析解的问题,蒙特·卡罗方法是一种有效的求出数值解的方法。一般蒙特·卡罗方法在数学中最常见的应用就是蒙特·卡罗积分1。
Monte Carlo Integration (蒙特卡罗积分)
- 蒙特卡罗积分法是使用蒙特卡罗方法求解定积分的一种应用
- 蒙特卡罗积分求解一般有两类方法: 2
- 投点法(频率法)
- 平均值法(期望法)
- 投点法
- 设要计算的定积分为
I
=
∫
a
b
g
(
x
)
d
x
I=\int_a^bg(x)dx
I=∫abg(x)dx, 其中
a
,
b
a,b
a,b为有限值, 被积函数
g
(
x
)
g(x)
g(x)是连续的随机变量
ξ
\xi
ξ的概率密度函数. 因此
g
(
x
)
g(x)
g(x)满足如下条件:
- g ( x ) g(x) g(x)非负
- ∫ − ∞ ∞ g ( x ) d x = 1 \int_{-\infty}^\infty g(x)dx = 1 ∫−∞∞g(x)dx=1
- 因此, I I I是一个概率积分, 其积分值等于概率 P ξ ( a ≤ ξ < b ) P_{\xi}(a\le \xi <b) Pξ(a≤ξ<b), 即 I = P ξ ( a ≤ ξ < b ) I = P_{\xi}(a\le \xi <b) I=Pξ(a≤ξ<b)
- 理解, 因为渲染方程求解的当前着色点在某个方向上的Radiance. 而在计算Radiance时, 能量总和一般不会大于1(能量守恒), 正好渲染方程符合蒙特卡罗积分性质, 其是一个概率密度函数(可能理解有误).
- 投点法求解: 给定分布
g
(
x
)
g(x)
g(x)的随机投点的办法, 求积分的实验步骤:
- 产生服从给定分布 g ( x ) g(x) g(x)的随机变量值 x ξ x_\xi xξ
- 检查 x ξ x_\xi xξ是否落入积分区域内, 即 a ≤ ξ < b a\le \xi <b a≤ξ<b时, 记录 x ξ x_\xi xξ落入积分区域一次
- 假设 N N N次实验后, x ξ x_\xi xξ落入积分区域的总次数为 m m m, 有: I ˉ = m N \bar{I} = \frac{m}{N} Iˉ=Nm作为该概率积分的近似值. 即, I ⋍ m N I \backsimeq \frac{m}{N} I⋍Nm
- 注意投点法求解的是概率积分. 最终的积分不会大于1. 感觉是在求 N N N次实验后随机变量落入积分区域的概率.
- 设要计算的定积分为
I
=
∫
a
b
g
(
x
)
d
x
I=\int_a^bg(x)dx
I=∫abg(x)dx, 其中
a
,
b
a,b
a,b为有限值, 被积函数
g
(
x
)
g(x)
g(x)是连续的随机变量
ξ
\xi
ξ的概率密度函数. 因此
g
(
x
)
g(x)
g(x)满足如下条件:
- 平均值法(期望法)
-
任取一组相互独立的、同分布的随机变量 { ξ i } \{\xi_i\} {ξi}, ξ i \xi_i ξi在 [ a , b ] [a,b] [a,b]中服从分布律 f ( x ) f(x) f(x), 令 g ∗ ( x ) = g ( x ) f ( x ) g^*(x) = \frac{g(x)}{f(x)} g∗(x)=f(x)g(x), 则 { g ∗ ( ξ i ) } \{g^*(\xi_i)\} {g∗(ξi)}也是一组相互独立、同分布的随机变量,而且:
E g ∗ ( ξ i ) = ∫ a b g ∗ ( x ) f ( x ) d x = ∫ a b g ( x ) d x = I (11) \tag{11} Eg^*(\xi_i) = \int_a^bg^*(x)f(x)dx = \int_a^bg(x)dx = I Eg∗(ξi)=∫abg∗(x)f(x)dx=∫abg(x)dx=I(11)
由大数定理
P r ( lim N → ∞ 1 N ∑ i = 1 N g ∗ ( ξ i ) = I ) = 1 (12) \tag{12}P_r(\lim_{N\to\infty}\frac{1}{N}\sum_{i=1}^Ng^*(\xi_i) = I) = 1 Pr(N→∞limN1i=1∑Ng∗(ξi)=I)=1(12)
若选
I = 1 N ∑ i = 1 N g ∗ ( ξ i ) (13) \tag{13}I = \frac{1}{N}\sum_{i=1}^Ng^*(\xi_i) I=N1i=1∑Ng∗(ξi)(13)
则 I ˉ \bar{I} Iˉ依概率1收敛于 I I I. 平均值法就是用 I ˉ \bar{I} Iˉ作为 I I I的近似值. -
假设要计算的积分有如下形式
I = ∫ a b g ( x ) d x (14) \tag{14} I = \int_a^bg(x)dx I=∫abg(x)dx(14)
其中被积函数 g ( x ) g(x) g(x)在区间 [ a , b ] [a,b] [a,b]内可积. 任意选择一个有简便方法可以进行抽样的概率密度函数 f ( x ) f(x) f(x), 使其满足下列条件:- f ( x ) ≠ 0 f(x) \not=0 f(x)=0 当 g ( x ) ≠ 0 g(x)\not=0 g(x)=0时 a ≤ x ≤ b a\le x \le b a≤x≤b
- ∫ a b f ( x ) d x = 1 \int_a^bf(x)dx = 1 ∫abf(x)dx=1
-
如果记
g ∗ ( x ) = { g ( x ) f ( x ) , f ( x ) ≠ 0 0 , f ( x ) = 0 (15) \tag{15} g^*(x) = \begin{cases} \frac{g(x)}{f(x)}, &f(x) \not=0 \\ 0 ,& f(x) = 0 \end{cases} g∗(x)={f(x)g(x),0,f(x)=0f(x)=0(15)
那么积分(14)可以写为
I = ∫ a b g ∗ ( x ) f ( x ) d x (16) \tag{16}I = \int_a^bg^*(x)f(x)dx I=∫abg∗(x)f(x)dx(16)
因而求积分的实验步骤是:- 产生服从分布律 f ( x ) f(x) f(x)的随机变量 x i ( i = 1 , 2 , 3... , N ) x_i(i=1, 2,3...,N) xi(i=1,2,3...,N)
- 计算均值
I
ˉ
=
1
N
∑
i
=
1
N
g
∗
(
x
i
)
\bar{I} = \frac{1}{N}\sum_{i=1}^{N}g^*(x_i)
Iˉ=N1∑i=1Ng∗(xi)
并用它作为 I I I的近似值, 即 I ≃ I ˉ I\simeq\bar{I} I≃Iˉ
-
如果 a 、 b a、b a、b为有限值,那么 f ( x ) f(x) f(x)可取作为均匀分布:
f ( x ) = { 1 b − a , a ≤ x ≤ b , 0 , 其 他 (17) \tag{17} f(x) = \begin{cases} \frac{1}{b-a}, & a \le x \le b , \\ 0, & 其他\end{cases} f(x)={b−a1,0,a≤x≤b,其他(17)
这时公式(14)变为
I = ( b − a ) ∫ a b g ( x ) 1 b − a d x (18) \tag{18} I = (b-a)\int_a^bg(x)\frac{1}{b-a}dx I=(b−a)∫abg(x)b−a1dx(18)
具体实验步骤如下:- 产生 [ a , b ] [a, b] [a,b]上的均匀分布随机变量 x i ( i = 1 , 2 , . . . , N ) x_i(i=1, 2, ..., N) xi(i=1,2,...,N)
- 计算均值 I ˉ = ( b − a ) N ∑ i = 1 N g ( x i ) \bar{I} = \frac{(b-a)}{N} \sum_{i=1}^Ng(x_i) Iˉ=N(b−a)∑i=1Ng(xi), 并用它作为 I I I的近似值, 即 I ≃ I ˉ I\simeq\bar{I} I≃Iˉ
-
利用蒙特卡罗积分解渲染方程
L
o
(
p
,
ω
o
)
=
L
e
(
p
,
ω
o
)
+
∫
Ω
+
L
i
(
p
,
ω
i
)
f
r
(
p
,
ω
i
,
ω
o
)
(
n
⋅
ω
i
)
d
ω
i
(6)
\tag6 L_o(p, \omega_o) =L_e(p, \omega_o) + \int_{\Omega^+}L_i(p, \omega_i)f_r(p, \omega_i, \omega_o)(n\cdot\omega_i)d\omega_i
Lo(p,ωo)=Le(p,ωo)+∫Ω+Li(p,ωi)fr(p,ωi,ωo)(n⋅ωi)dωi(6)
渲染方程中去掉自发光项,即反射方程如下:
L
o
(
p
,
ω
o
)
=
∫
Ω
+
L
i
(
p
,
ω
i
)
f
r
(
p
,
ω
i
,
ω
o
)
(
n
⋅
ω
i
)
d
ω
i
(19)
\tag{19} L_o(p, \omega_o) = \int_{\Omega^+}L_i(p, \omega_i)f_r(p, \omega_i, \omega_o)(n\cdot\omega_i)d\omega_i
Lo(p,ωo)=∫Ω+Li(p,ωi)fr(p,ωi,ωo)(n⋅ωi)dωi(19)
- 从渲染方程(公式19)可以看出。该积分为公式(14)的形式, 其中 g ( x ) = L i ( p , ω i ) f r ( p , ω i , ω o ) ( n ⋅ ω i ) g(x) = L_i(p, \omega_i)f_r(p, \omega_i, \omega_o)(n\cdot\omega_i) g(x)=Li(p,ωi)fr(p,ωi,ωo)(n⋅ωi)
- 如果要求解公式(19), 可以根据公式(16)的求解过程将积分改写为:
I ˉ = 1 N ∑ i = 1 N g ∗ ( x i ) = 1 N ∑ i = 1 N g ( x i ) f ( x i ) = 1 N ∑ i = 1 N L i ( p , ω i ) f r ( p , ω i , ω o ) ( n ⋅ ω i ) f ( ω i ) (20) \tag{20} \bar{I} = \frac{1}{N}\sum_{i=1}^{N}g^*(x_i) = \frac{1}{N}\sum_{i=1}^{N}\frac{g(x_i)}{f(x_i)} = \frac{1}{N}\sum_{i=1}^{N}\frac{L_i(p, \omega_i)f_r(p, \omega_i, \omega_o)(n\cdot\omega_i)}{f( \omega_i)} Iˉ=N1i=1∑Ng∗(xi)=N1i=1∑Nf(xi)g(xi)=N1i=1∑Nf(ωi)Li(p,ωi)fr(p,ωi,ωo)(n⋅ωi)(20)
由于我们求解的是着色点上半球接收到所有Radiance, 因此概率密度函数(PDF)可以设为 f ( x ) = 1 / 2 π f(x) = 1/2\pi f(x)=1/2π(半球为 2 π 2\pi 2π, 均匀采样). 因此公式(20)改写为:
I ˉ = 1 N ∑ i = 1 N g ∗ ( x i ) = 1 N ∑ i = 1 N g ( x i ) f ( x i ) = 1 N ∑ i = 1 N L i ( p , ω i ) f r ( p , ω i , ω o ) ( n ⋅ ω i ) 1 / 2 π (21) \tag{21} \bar{I} = \frac{1}{N}\sum_{i=1}^{N}g^*(x_i) = \frac{1}{N}\sum_{i=1}^{N}\frac{g(x_i)}{f(x_i)} = \frac{1}{N}\sum_{i=1}^{N}\frac{L_i(p, \omega_i)f_r(p, \omega_i, \omega_o)(n\cdot\omega_i)}{1/2\pi} Iˉ=N1i=1∑Ng∗(xi)=N1i=1∑Nf(xi)g(xi)=N1i=1∑N1/2πLi(p,ωi)fr(p,ωi,ωo)(n⋅ωi)(21)
(在求解渲染方程代码中, PDF使用的是 1 / 4 π r 2 1/4\pi r^2 1/4πr2, 感觉更合理些. 但是PDF的选择是多样的, 只需满足 ∫ a b f ( x ) d x = 1 \int_a^bf(x)dx = 1 ∫abf(x)dx=1 并且在积分区间内 f ( x ) ≠ 0 f(x)\not=0 f(x)=0即可 ) - 最后根据路径追踪(当 n = 1 n=1 n=1时, 即路径追踪), 递归求解渲染方程即可得到当前着色点的最终的出射Radiance(如果自发光物体需要加上自发光项).
- 基础内容:递归的终止条件、光线的产生和以及直接光照和间接光照区分(略)
小结
复习渲染方程相关的理论知识, 更加深入的去理解渲染方程背后的数学原理., 温故而知新.