大気散乱 その1

大気散乱

フォグを実装したときに, 大気散乱としてレイリー散乱とミー散乱の位相関数を利用して,フォグの濃度を計算しました。 これをもう少しちゃんと考えてみて,スカイボックスを自動生成してみたいと思います。
"その1"と付けたのは,まだ実装とか調査途中だからです。"その2"があるかは分かりません。

http://amd-dev.wpengine.netdna-cdn.com/wordpress/media/2012/10/ATI-LightScattering.pdf
https://ebruneton.github.io/precomputed_atmospheric_scattering/
https://www.guerrilla-games.com/read/decima-engine-advances-in-lighting-and-aa
http://www.crytek.com/cryengine/presentations/real-time-atmospheric-effects-in-games-revisited

散乱

光が大気中の粒子と,周囲の光源の散乱によってどの様に光が見えるかを考えます。 ATIのドキュメントや,Precomputed Atmospheric Scatteringのペーパーを読んで考えてみます。

f:id:hikita12312:20180209000343p:plain:w300

{p_0}の観測者が{p_1}から発せられる光{L_{\mathrm{bg}}}を見ます。 この時,観測者に届く光{L(p_0,p1)}は,

{ \displaystyle
\begin{eqnarray}
L({\bf p}_0, {\bf p}_1) &=& L_{\mathrm{bg}} f_{\mathrm{ex}}({\bf p}_0, {\bf p}_1) + \int d{\bf p} f_{\mathrm{ex}}({\bf p}_0, {\bf p}) \sum_{R,M} \beta_{R,M}({\bf p})\int_{\Omega} d\omega F_{R,M}(\theta) L_{R,M}(\theta,\phi)
\end{eqnarray}
}

{f_{\mathrm{ex}}}は光が届くまでに減衰する成分です。{R,M}はそれぞれレイリー散乱とミー散乱による項を示します。 {\beta(p)}は地点{p}における粒子の濃度によるスケール成分で,{F(\theta)}は光源と観測者方向のベクトルのなす角{\theta}によるレイリー散乱とミー散乱の位相関数です。

{ \displaystyle
\begin{eqnarray}
F_{\mathrm {Reyligh}}(\theta) &=& \frac{3}{16 \pi}(1+\cos^2{\theta}) \\
F_{\mathrm {Mie}}(\theta) &=& \frac{1}{4 \pi}\frac{(1+g)^2}{(1+g^2-2g\cos{\theta})^{3/2}}
\end{eqnarray}
}

結局のところ,

  • もともとの背景色が粒子によって減衰して届く光
  • 太陽光が各位置の粒子によって散乱されて届く光

の加算されたものが,目に見える光となるでしょう。さらにここから散乱された光がさらに散乱される多重散乱も考えることが出来ますが,簡単のため無視します。

さて,グーグルで調べると,地球から太陽までの距離は149,600,000 kmですが,粒子が存在する大気の厚さは500km程度だそうです。 散乱される場所から,太陽方向と観測者方向のなす角は,積分範囲では変化は微小と考えても良いでしょう。 また,多重散乱も考えないので,太陽以外に光源が無いとして,さらに太陽項の存在範囲は立体角では微小な領域{\Delta\omega}で一様な光であるとみなします。 {f_{\mathrm{ex}}(p_0, p)}も散乱光に対してはそれほど減衰の効果が無いとして,1に近似してしまいます。(※ これは本当はダメな気がする...)

{ \displaystyle
\begin{eqnarray}
L({\bf p}_0, {\bf p}_1) &=& L_{\mathrm{bg}} f_{\mathrm{ex}}({\bf p}_0, {\bf p}_1) + \sum_{R,M} F_{R,M}(\theta) L_{R,M} \Delta\omega \int d{\bf p} \beta_{R,M}({\bf p})
\end{eqnarray}
}

ここで{\int dp \beta_{R,M}(p)}積分は,{p_0}から{p_1}までの光路で吸収されなかった光が散乱された割合と考えても良い気がしてきます。 {f_{\mathrm{ex}}}を光の透過率と考えて {f_{\mathrm{ex}}=T}として,

{ \displaystyle
\begin{eqnarray}
L({\bf p}_0, {\bf p}_1) &=& L_{\mathrm{bg}} T + \sum_{R,M} F_{R,M}(\theta) L_{R,M} \Delta\omega (1-T_{R,M})
\end{eqnarray}
}

つまり,以前実装したフォグの実装は, フォグ濃度が{1-T}に対応していて,背景色とフォグカラーを適度に補間するという考えは,あながち間違いではなかったのかもしれません。

球面上のフォグ透過率の積分

フォグ濃度の計算をしたときと同様に,高さに対して指数的に粒子の濃度が変わっていく高さフォグを仮定します。

光が単位距離だけ進んだら,{a(y)=a_0 \exp{(-b y)}}倍だけ減衰するときに, {(0,y_0)}の位置から,{(d,y_0+h)}の位置までの積分結果は以下のようになるのでした。

{ \displaystyle
\begin{eqnarray}
T &=& \exp{\Big\{ \frac{a_0 d \exp{(-b y_0)} }{b h} \Big( \exp{(-b h)} - 1\Big) \Big\} }
\end{eqnarray}
}

今回は地球の球面上での高さフォグとなるので,単純には経路積分できません。
とりあえず,経路を分割して,各領域では球面ではなく,平坦で上の高さフォグの積分計算が利用できると仮定して積分します。
(もしかしたら厳密解があるのかもしれませんが...)

f:id:hikita12312:20180210103357p:plain:w300
/// @brief 球面上の高さフォグの積分計算
/// @param h+R 観測者の高さ + 地球の半径
/// @param integralViewVec 視点方向ベクトル,長さが積分範囲
/// @param R 地球の半径
/// @brief a0 高さフォグの地表における濃度 (reyleigh,mie)
/// @brief b 高さフォグの単位高さあたりの減衰率 (reyleigh,mie)
float2 GetTransmittance(float hR,float3 integralViewVec,float2 a0,float2 b,float R,uint N)
{
    // 1stepごとの位置差分
    const float3 step = integralViewVec / (N - 1);
    const float stepLength2 = dot(step, step);
    // 地球中心から観測者までのベクトル
    const float3 P0 = float3(0.0f, hR, 0.0f);
    //// 積分計算
    float2 sum = 0.0f;
    [loop]
    for (uint i = 0; i < N; ++i)
    {
        // 積分の区分範囲を計算
        const float3 Q1 = P0 + step * i;
        const float3 Q2 = P0 + step * (i + 1);
        // Q1,Q2の地表からの高さを計算
        const float h1 = max(0.0f, length(Q1) - R);
        const float h2 = max(0.0f, length(Q2) - R);
        // 球面座標系のQ1,Q2の垂直/水平距離
        const float deltaH = (h2 - h1);
        const float deltaT = sqrt(max(stepLength2 - deltaH*deltaH, 0.0f));
        // 球面座標系における高さFogの経路積分を直線とみなして計算
        const float invK = deltaT / (deltaH + 0.0001f); // 傾きの逆数(発散封じの項付き)
        sum += exp(-b*h1) * (exp(-b*deltaH) - 1.0f) * invK;
    }
    return exp(sum * a0 / b);
}

この透過率を使って,大気散乱を計算します。

結果

256x256サイズのキューブマップテクスチャと対応する視線方向の角度の大気散乱の色を計算して,スカイボックス上に投影した結果です。 計算的に地表となる地点には適当に地面っぽい色を置いています。 また,太陽は256x256のテクスチャにレンダリングすると荒すぎるので別途レンダリングしています。

f:id:hikita12312:20180210103417p:plain:w600
f:id:hikita12312:20180210103438p:plain:w600

大気っぽい気はしますが,微妙な感じはします。
パラメータ調整が甘いんでしょうか。多重散乱なども真面目に計算するべきなのでしょうか。 そのうち直したいと思います。