J1N9H3
雾、云、天空和大气的环境光源

天空和大气的渲染 单散射实现

参考 1993 年工作的大气渲染单次散射实现

完整实现代码

Alpha3D v0.1.6

Image 1

参照论文 Nishita 1993 - Display of The Earth Taking into Account Atmospheric Scattering 实现。

天空和大气的渲染

观察射线从观察点 PvP_v 出发,在 PaP_a 处进入大气层,并在 PbP_b 处与地球表面相交。摄像机最终接收到的光,可以描述为路径 PaPbP_aP_b 上各点产生的大气散射光 IvI_v 的累积。太阳光到达路径上的某个大气采样点 PP 后,一部分光被散射到观察方向,并在从 PP 传播到观察点 PvP_v 的过程中继续受到大气消光作用。

地球大气主要由氮气、氧气等气体分子,以及臭氧和气溶胶等成分组成。臭氧等物质会选择性地吸收特定波长的光;气体分子和气溶胶则会使光发生散射。气体分子的尺度远小于可见光波长,其散射通常使用瑞利散射描述;气溶胶颗粒的尺度与可见光波长接近,其散射通常使用米氏散射描述。

光的散射由散射系数和相位函数共同描述,光的吸收与整体衰减则由吸收系数以及比尔–朗伯定律描述。摄像机接收到的大气散射光,是观察路径上各采样点产生的散射贡献经过两段路径衰减后累积得到的结果。

光散射的描述

定义光的散射模型:输入波长为 λ\lambda 的入射光强度,输出 θ\theta 角度散射后的光强度。

I(λ,θ)=I0(λ)KρFR(θ)λ4.(1)I(\lambda,\theta)=I_0(\lambda)\,K\rho\,\frac{F_R(\theta)}{\lambda^4}. \tag{1}

其由 KKρ\rhoFR(θ)F_R(\theta) 描述。

  • KK 与空气分子物理性质有关,为常数;
  • ρ\rho 描述当前高度的相对空气密度,ρ(h)=exp(hH0)\rho(h)=\exp\left(-\frac{h}{H_0}\right)
  • FR(θ)F_R(\theta) 表示散射相函数,描述散射的方向分布。

指数密度模型在 Shader 中可以直接写成:

其中返回值的两个分量分别保存瑞利粒子和米氏粒子的相对密度。后续只需乘以各自的散射系数,就能得到采样点处的散射强度。

光学深度

定义光沿整条路径传播时,沿途每一小段大气造成的消光作用之和。

考虑路径上的一小段 dsds。该处的大气密度比例是 ρ(s)\rho(s),标准大气下单位长度的消光系数是 β(λ)\beta(\lambda),这一小段产生的光学深度为 dt=β(λ)ρ(s)dsdt=\beta(\lambda)\rho(s)\,ds

t(S,λ)=0Sβ(λ)ρ(s)ds=4πKλ40Sρ(s)ds.(2)t(S,\lambda) = \int_0^S \beta(\lambda)\rho(s)\,ds = \frac{4\pi K}{\lambda^4} \int_0^S\rho(s)\,ds. \tag{2}

连续积分在 GPU 上需要离散为有限个采样点。将长度为 rayLength 的路径均分为 NN 段,并使用中点法累积柱密度:

columnDensity 对应 ρ(s)ds\int\rho(s)\,ds,最后分别乘以 RGB 消光系数,得到式(2)中的光学深度 τ\tau

Beer–Lambert 衰减

光穿过介质时,光强会随着传播距离按指数形式减小

Iout=Iineτ.(3)I_{\mathrm{out}}=I_{\mathrm{in}}e^{-\tau}. \tag{3}

其中 τ\tau 表示光学深度,eτe^{-\tau} 即透射率(transmittance),表示剩余的光比例。

GLSL 的 exp 会逐分量计算,因此 RGB 三个波段可以同时应用 Beer–Lambert 衰减。

单次散射的路径积分

假设 Is(λ)I_s(\lambda) 为平行光,因此 θ\theta 为常数,描述 PP 点的光照强度 Ip(λ)I_p(\lambda) 如下:

Ip(λ)=Is(λ)KρFR(θ)λ4exp[t(PPc,λ)].(4)I_p(\lambda)=I_s(\lambda)\,K \rho\,\frac{F_R(\theta)}{\lambda^4}\exp\left[-t(PP_c,\lambda)\right]. \tag{4}

其中,光学深度描述如下:

t(PPc,λ)=PcPβ(l,λ)ρ(l)dl.(5)t(PP_c,\lambda)=\int_{P_c}^{P}\beta(l,\lambda)\rho(l)\,dl. \tag{5}

由 Beer–Lambert 可以得到 PP 点的光照强度:

Ipv(λ)=Ip(λ)exp[t(PPa,λ)].(6)I_{pv}(\lambda)=I_p(\lambda)\exp\left[-t(PP_a,\lambda)\right]. \tag{6}

将式(4)代入式(6),可以得到采样点 PP 对观察点 PvP_v 的光强贡献:

Ipv(λ)=Is(λ)Kρ(P)FR(θ)λ4exp[t(PPc,λ)t(PPa,λ)].(7)I_{pv}(\lambda)=I_s(\lambda)\,K\rho(P)\,\frac{F_R(\theta)}{\lambda^4}\exp\left[-t(PP_c,\lambda)-t(PP_a,\lambda)\right]. \tag{7}

累计路径 PbPaP_bP_a 上的每个采样点的光强贡献,进行积分:

Iv(λ)=PbPaIpv(P,λ)ds.(8)I_v(\lambda)=\int_{P_b}^{P_a}I_{pv}(P,\lambda)\,ds. \tag{8}

将式(7)代入式(8),得到:

Iv(λ)=Is(λ)KFR(θ)λ4PbPaρ(P)exp[t(PPc,λ)t(PPa,λ)]ds.(9)I_v(\lambda)=I_s(\lambda)\frac{K F_R(\theta)}{\lambda^4}\int_{P_b}^{P_a}\rho(P)\exp\left[-t(PP_c,\lambda)-t(PP_a,\lambda)\right]\,ds. \tag{9}

根据散射系数的定义:

βR(λ)=4πKλ4,(10)\beta_R(\lambda)=\frac{4\pi K}{\lambda^4}, \tag{10}

再定义归一化相函数:

pR(θ)=FR(θ)4π,(11)p_R(\theta)=\frac{F_R(\theta)}{4\pi}, \tag{11}

则有:

βR(λ)pR(θ)=KFR(θ)λ4.(12)\beta_R(\lambda)p_R(\theta)=\frac{K F_R(\theta)}{\lambda^4}. \tag{12}

因此,式(9)可以改写为:

Iv(λ)=Is(λ)βR(λ)pR(θ)PbPaρR(P)exp[t(PPc,λ)t(PPa,λ)]ds.(13)I_v(\lambda)=I_s(\lambda)\beta_R(\lambda)p_R(\theta)\int_{P_b}^{P_a}\rho_R(P)\exp\left[-t(PP_c,\lambda)-t(PP_a,\lambda)\right]\,ds. \tag{13}

代入光学深度公式 (5),得到:

Iv(λ)=Is(λ)βR(λ)pR(θ)PbPaρR(P)exp[PPcβ(l1,λ)ρ(l1)dl1PPaβ(l2,λ)ρ(l2)dl2]ds.(14)I_v(\lambda) = I_s(\lambda)\beta_R(\lambda)p_R(\theta) \int_{P_b}^{P_a} \rho_R(P) \exp\left[ -\int_{P}^{P_c} \beta(l_1,\lambda)\rho(l_1)\,\mathrm{d}l_1 -\int_{P}^{P_a} \beta(l_2,\lambda)\rho(l_2)\,\mathrm{d}l_2 \right] \,\mathrm{d}s. \tag{14}

式(14)对应两层积分:外层沿观察射线寻找产生散射的采样点,内层计算太阳到采样点以及采样点到观察者两段路径的透射率。离散后的核心结构如下:

其中 sunTransmittance 对应 exp[t(PPc)]\exp[-t(PP_c)]viewTransmittance 对应 exp[t(PPa)]\exp[-t(PP_a)]localScattering * stepLength 则是采样点 PP 对外层积分的一小段贡献。

大气密度分布

实现分别计算瑞利散射粒子、米氏散射粒子和吸收介质的相对密度。设采样点海拔高度为 hh,瑞利和米氏密度都使用指数大气模型:

ρR(h)=exp(hHR),ρM(h)=exp(hHM).(12)\rho_R(h)=\exp\left(-\frac{h}{H_R}\right), \qquad \rho_M(h)=\exp\left(-\frac{h}{H_M}\right). \tag{12}

HRH_RHMH_M 分别是瑞利与米氏密度的标高。标高越大,密度随高度下降得越慢。气溶胶通常集中在近地面,因此米氏标高通常明显小于瑞利标高。

吸收介质使用以 hAh_A 为中心、由 HAH_A 控制宽度的钟形分布:

o(h)=hAhHA,ρA(h)=ρR(h)o(h)2+1.(13)o(h)=\frac{h_A-h}{H_A}, \qquad \rho_A(h)=\frac{\rho_R(h)}{o(h)^2+1}. \tag{13}

分母让密度在吸收层中心附近较大,并随偏离中心的距离平滑下降;额外乘以瑞利密度,使高空部分继续随整体大气密度衰减。这是适合实时渲染的简化分布,而不是对真实臭氧浓度剖面的严格拟合。

可见光光谱

可见光的波长大约位于 380780nm380\text{–}780 \, \mathrm{nm}。下图描述了波长从短到长,颜色的变化:

天空和大气的渲染

整体关系可以表示为:

短波长紫、蓝、绿、黄、橙、红长波长\text{短波长} \rightarrow \text{紫、蓝、绿、黄、橙、红} \rightarrow \text{长波长}

瑞利散射

瑞利散射(Rayleigh)描述空气分子对光的散射,主要决定天空、大气边缘以及远距离景物的基础颜色。其散射强度近似与波长的四次方成反比:

βR(λ)1λ4\beta_R(\lambda)\propto\frac{1}{\lambda^4}

因此,波长较短的蓝光比波长较长的红光更容易被散射。

白天,天空蓝光被更多地散射到观察者方向,使天空呈蓝色。

NoneNone
RayleighRayleigh

日出和日落时,太阳光经过更长的大气路径,蓝光被大量散射出去,直射太阳光偏橙红色。

NoneNone
RayleighRayleigh

瑞利散射相函数通常写为:

pR(θ)=316π(1+cos2θ)p_R(\theta)=\frac{3}{16\pi}\left(1+\cos^2\theta\right)

它表示散射光在前向和后向较强,在侧向相对较弱,并且前后方向基本对称。

实现时直接传入 cosθ\cos\theta,避免在每个采样点调用 acos。由于太阳光被视为平行光,同一条观察射线上的 cosTheta 保持不变,因此相函数可以移到积分循环外计算一次。

米氏散射

米氏散射(Mie)用于描述气溶胶、灰尘、水滴等较大颗粒对光的散射。米氏散射通常具有明显的前向散射特征,即光更容易沿原来的传播方向被散射。因此,朝太阳方向观察时会看到更强的亮度和光晕。

通常使用 Henyey–Greenstein 相函数近似其方向分布:

pM(θ)=1g24π(1+g22gcosθ)3/2p_M(\theta) = \frac{1-g^2} {4\pi\left(1+g^2-2g\cos\theta\right)^{3/2}}

max 用来防止 ggcosθ\cos\theta 接近 1 时分母过小。这个方向正是米氏前向散射最集中的区域,也对应太阳周围明亮的光晕。

其中 gg 是米氏各向异性系数:

gg散射方向
=0= 0各向同性散射
>0> 0前向散射
<0< 0后向散射
1\rightarrow 1高度集中的前向散射

实现中 g=0.7g=0.7,表示具有比较明显的前向散射,所以太阳周围会形成较集中的亮光晕。

RayleighRayleigh
Rayleigh + MieRayleigh + Mie

三个颜色通道受到近似相同强度的散射,Mie 散射通常表现为白色或灰白色,而不像 Rayleigh 散射那样明显偏蓝。

臭氧吸收

臭氧层选择性地削弱部分可见光,臭氧对绿光吸收最强,对红光次之,对蓝光最弱(βG>βR>βB\beta_G \gt \beta_R \gt \beta_B)。

实现使用一个以 absorptionHeight 为中心的平滑密度层近似臭氧分布:

计算消光时,臭氧密度只乘以吸收系数,不会像瑞利或米氏粒子一样产生朝观察方向的散射贡献:

日出、日落时,光在大气中的传播路径变长,臭氧对绿色和黄绿色成分的选择性吸收更加明显,从而改变天空的色彩平衡,使部分区域呈现更明显的蓝紫色或粉紫色过渡。

Rayleigh + MieRayleigh + Mie
Rayleigh + Mie + AbsorptionRayleigh + Mie + Absorption

从太空观察时,大气分子产生的瑞利散射使行星边缘呈现蓝色光环,臭氧吸收会进一步削弱部分红光和绿光,从而改变并强化大气边缘的颜色表现。

Rayleigh + MieRayleigh + Mie
Rayleigh + Mie + AbsorptionRayleigh + Mie + Absorption

On this page