大气散射（七） 大气散射的shader

Mr_Swilder
457 阅读5分钟

一、采样视线方向

让我们回顾一下我们最近推导出的大气散射方程：

I=IS∑P∈AB‾S(λ,θ,h)T(CP‾)T(PA‾)dsI= I_S \sum_{P \in \overline{AB}} {S\left(\lambda, \theta, h\right) T\left(\overline{CP}\right) T\left(\overline{PA}\right) ds }

我们接收到的光量等于来自太阳的光量 ISI_S，乘以 AB‾\overline{AB} 中每个点 P 的单独贡献的总和。

我们可以直接在我们的着色器中实现这个函数。然而，有一些优化可以做。在之前的教程中，曾经暗示过可以进一步简化这个表达式。我们可以将散射函数分解为其两个基本组成部分：

S(λ,θ,h)=β(λ,h)γ(θ)=β(λ)ρ(h)γ(θ)S \left(\lambda, \theta, h\right ) = \beta \left(\lambda, h \right ) \gamma\left(\theta\right) = \beta \left(\lambda\right )\rho\left(h\right) \gamma\left(\theta\right)

相位函数 γ(θ)\gamma\left(\theta\right) 和海平面处的散射系数 β(λ)\beta \left(\lambda\right ) 是与求和无关的常数，因为角度 θ\theta 和波长λ\lambda 不依赖于采样点。因此，它们可以被提取出来：

I=IS β(λ)γ(θ)∑P∈AB‾T(CP‾)T(PA‾)ρ(h)dsI = I_S \, \beta \left(\lambda\right ) \gamma\left(\theta\right) \sum_{P \in \overline{AB}} { T\left(\overline{CP}\right) T\left(\overline{PA}\right) \rho\left(h\right) ds }

这个新表达式在数学上等价于以前的表达式，但计算效率更高，因为一些最重要的部分已经被提取出来。

我们现在可以开始实现它了。我们应该考虑无限多个点 P。对 I 的一个合理近似是将 AB‾\overline{AB} 分成几个长度为 ds 的较小段，并累积每个单独段的贡献。这样做的时候，我们假设每个段都足够小，以至于其密度是恒定的。一般来说，这并不是真实情况，但如果 ds 足够小，我们仍然可以得到一个相当好的近似。

image.png AB‾\overline{AB} 中的段数称为视线采样点，因为所有段都位于视线射线上。在着色器中，这将是 _ViewSamples 属性。通过将其作为属性，它可以从材质检视器中访问。这允许我们在性能方面减少着色器的精度。

下面的代码段允许循环遍历大气中的所有段。

// 数值积分以计算
// AB 中每个点 P 的光贡献
float3 totalViewSamples = 0;
float time = tA;
float ds = (tB-tA) / (float)(_ViewSamples);
for (int i = 0; i < _ViewSamples; i ++)
{
    // 点的位置
    //（在视线采样段的中间采样）
    float3 P = O + D * (time + ds * 0.5);
    // T(CP) * T(PA) * ρ(h) * ds
    totalViewSamples += viewSampling(P, ds);
    time += ds;
}
// I = I_S * β(λ) * γ(θ) * totalViewSamples
float3 I = _SunIntensity *  _ScatteringCoefficient * phase * totalViewSamples;

变量 time 用于跟踪我们离起始点 O 有多远，并在每次迭代后增加 ds。

二、光学深度PA

沿着视线 (\overline{AB}) 的每个点都对我们绘制的像素的最终颜色贡献了自己的部分。从数学上讲，这个贡献是求和符号中的数量：

I=IS β(λ)γ(θ)∑P∈AB‾T(CP‾)T(PA‾)ρ(h)ds⏟light contribution of L(P)I = I_S \, \beta \left(\lambda\right ) \gamma\left(\theta\right) \sum_{P \in \overline{AB}} \underset{\text{light contribution of}\,L\left(P\right)}{\underbrace{T\left(\overline{CP}\right) T\left(\overline{PA}\right) \rho\left(h\right) ds}}

像在上一段中所做的那样，让我们尝试进一步简化它。我们可以通过用其实际定义替换 T 来进一步扩展上述表达式：

T(XY‾)=exp⁡{−β(λ)D(XY‾)}T\left(\overline{XY}\right) =\exp\left\{ - \beta\left(\lambda\right) D\left(\overline{XY}\right) \right\}

CP‾\overline{CP}和 PA‾\overline{PA} 上的透射率的乘积变成了：

T(CP‾)T(PA‾)=T\left(\overline{CP}\right) T\left(\overline{PA}\right)=

=exp⁡{−β(λ)D(CP‾)}⏟T(CP‾) exp⁡{−β(λ)D(PA‾)}⏟T(PA‾)==\underset{T\left(\overline{CP}\right) }{\underbrace{ \exp\left\{- \beta\left(\lambda\right) D\left(\overline{CP}\right)\right \} }} \, \underset{T\left(\overline{PA}\right) }{\underbrace{ \exp\left\{- \beta\left(\lambda\right) D\left(\overline{PA}\right) \right \} }}=

=exp⁡{−β(λ)(D(CP‾)+D(PA‾))}= \exp\left\{- \beta\left(\lambda\right) \left( D\left(\overline{CP}\right) + D\left(\overline{PA}\right) \right) \right \}

联合透射率被建模为指数衰减，其系数是光线（CP‾\overline{CP} 和 PA‾\overline{PA}）的路径上的光学深度之和，乘以海平面处的散射系数（β\beta，其中 h=0h=0）。

我们首先要计算的量是段 PA‾\overline{PA} 的光学深度，它从进入大气的点穿过大气，直到我们在 for 循环中当前正在采样的点。让我们回顾一下光学深度的定义：

D(PA‾)=∑Q∈PA‾exp⁡{−hQH} dsD\left( \overline{PA}\right)=\sum_{Q \in \overline{PA}} { \exp\left\{-\frac{h_Q}{H}\right\} } \, ds

如果要天真地实现这一点，我们将在一个循环中对 P 和 A 之间的点进行采样。这是可能的，但效率非常低下。实际上，D(PA‾)D\left( \overline{PA}\right) 就是我们已经分析的外层 for 循环中当前段的光学深度。如果我们计算当前以 P 为中心的段的光学深度（opticalDepthSegment），并在 for 循环中持续累积它（opticalDepthPA），我们可以节省大量计算。

// 光学深度累加器
float opticalDepthPA = 0;
// 数值积分以计算
// AB 中每个点 P 的光贡献
float time = tA;
float ds = (tB-tA) / (float)(_ViewSamples);
for (int i = 0; i < _ViewSamples; i ++)
{
    // 点的位置
    //（在视线采样段的中间采样）
    float3 P = O + D * (time + viewSampleSize*0.5);
    // 当前段的光学深度
    // ρ(h) * ds
    float height = distance(C, P) - _PlanetRadius;
    float opticalDepthSegment = exp(-height / _ScaleHeight) * ds;
    // 累加
    // 累加光学深度 opticalDepthPA += opticalDepthSegment; ...  
    time += ds;
}

三、光线采样

如果我们回顾一下点P的光贡献的表达式，我们会发现唯一需要的量是线段(\overline{CP})的光学深度：

L(P)=exp⁡{−β(λ)(D(CP‾)+D(PA‾))}⏟合并透射率 ρ(h)ds⏟线段的光学深度L\left(P\right) = \underset{\text{合并透射率}}{\underbrace{\exp\left\{- \beta\left(\lambda\right) \left( D\left(\overline{CP}\right) + D\left(\overline{PA}\right) \right) \right \}}} \, \underset{\text{线段的光学深度}}{\underbrace{\rho\left(h\right) ds}}

我们将计算线段(\overline{CP})的光学深度的代码移到一个名为lightSampling的函数中。这个名称来自于光线，它是从P开始并指向太阳的线段。我们称它穿出大气层的点为C。

然而，lightSampling函数不仅仅计算(\overline{CP})的光学深度。到目前为止，我们只考虑了大气层的贡献，忽略了实际行星的作用。我们的方程没有考虑到从P向太阳发出的光线可能会击中行星。如果这种情况发生，到目前为止所做的所有计算都必须被丢弃，因为没有光实际上会到达相机。

image.png 在上面的图中，很容易看出应该忽略P0P_0的光贡献，因为太阳光没有照到P0P_0。在循环遍历从P到C之间的点时，lightSampling函数还会检查行星是否被击中。这可以通过检查点的高度是否为负来完成。

bool lightSampling
(    float3 P,    // Current point within the atmospheric sphere
    float3 S,    // Direction towards the sun
    out float opticalDepthCA
)
{
    float _; // don't care about this one
    float C;
    rayInstersect(P, S, _PlanetCentre, _AtmosphereRadius, _, C);

    // Samples on the segment PC
    float time = 0;
    float ds = distance(P, P + S * C) / (float)(_LightSamples);
    for (int i = 0; i < _LightSamples; i ++)
    {
        float3 Q = P + S * (time + lightSampleSize*0.5);
        float height = distance(_PlanetCentre, Q) - _PlanetRadius;
        // Inside the planet
        if (height < 0)
            return false;

        // Optical depth for the light ray
        opticalDepthCA += exp(-height / _RayScaleHeight) * ds;

        time += ds;
    }

    return true;
}

该函数首先使用rayIntersect计算点C。然后，它将线段PA‾\overline{PA}分成长度为ds的_LightSamples个片段。光学深度的计算与最外层循环中使用的计算相同。

如果行星被击中，该函数返回false。我们可以使用这一点来更新最外层循环中缺失的代码，将“....”替换为以下代码：

// D(CP)
float opticalDepthCP = 0;
bool overground = lightSampling(P, S);
if (overground)
{
    // Combined transmittance
    // T(CP) * T(PA) = T(CPA) = exp{ -β(λ) [D(CP) + D(PA)]}
    float transmittance = exp
    (
        -_ScatteringCoefficient *
        (opticalDepthCP + opticalDepthPA)
    );
    // Light contribution
    // T(CPA) * ρ(h) * ds
    totalViewSamples += transmittance * opticalDepthSegment;
}

现在，我们已经考虑了所有要素，我们的着色器已经完成。

avatar
文章
阅读
粉丝