A Radiative Transfer Framework for Spatially-Correlated Materials

[1805.02651] A Radiative Transfer Framework for Spatially-Correlated Materials

f:id:mizuooon:20180605162830p:plain

SIGGRAPH 2018 の論文。

ボリュームレンダリングに使用する放射輸送方程式 (RTE) は媒質内の粒子が空間的に相関性のない配置をしていることを前提としていたが、 実はこの粒子の空間的な相関というのはボリュームの見た目に大きな影響を及ぼす。 この論文では、 粒子間の空間的相関性を考慮した RTE は Generalized Boltzmann Equation (GBE) として他分野で提案されていたものを レンダリングフレームワーク上で適用できるように拡張し、 粒子間の空間的相関性を考慮したボリュームレンダリングを可能にすることを目的としている。

ランベルト・ベールの法則

まずは、粒子の空間的相関性を考慮しない従来のモデルを考える。 光がある一様な媒質の中を直進する際、消散係数(extinction coefficient)を定数  \mu とすると その減衰は

 \displaystyle \frac{dL(\boldsymbol{x})}{ds} = -\mu L(\boldsymbol{x})

として表される。 距離 t を進むうち輝度 L の減衰量は、これを解いて

 (1 - exp(-\mu t)) L(\boldsymbol{x})

となり、つまり  1 - exp(-\mu t) の割合の光が失われ、 -exp(-\mu t) の割合の光が透過する。 (ちなみに、この光が透過した場合の式は一度も媒質中の粒子に衝突しなかったということであり、 k=0 の特殊なポアソン過程とも考えられる。)

このように光の減衰が距離に関して指数的 (ポアソン的) に起こるのをランベルト・ベールの法則と呼ぶ。

粒子の空間的相関性

以上が粒子に空間的相関性が無い場合についてであるが、相関性を考えると上記のモデルが正しくなくなる。

まずは空間的相関性が何かについて述べる。 これは単純に言うとある粒子の位置が別の粒子の位置に対して情報量を持っているか、を表す。 相関性が無い場合は、個別の粒子はそれぞれ完全にランダムに配置される(上図左)。 一方、相関性が正の場合は粒子の配置には空間的な偏りが生じ(上図中央)、 また相関性が負の場合はあるパターンの上で規則的に配置されることになる(上図右)。 このとき、 正の相関ではある粒子の近くには別の粒子があることが推測でき、 同様に負の相関でもパターンを利用して別の粒子の位置が推測できる。 つまり粒子の位置が別の粒子の位置に対する情報量を持っていることがわかる。 負の相関が規則正しい配列を示す、というのはやや直感に反するが、 ある粒子の近傍には他の粒子が無いのを表している。多分。

粒子間の空間的相関性の影響

こうした空間的相関性を持った粒子を考えると、 正相関がある場合は粒子同士の重なり合いが大きくなるため、 実効的な消散効果が指数的な減衰よりも小さくなる。 その結果、ボリュームは無相関の媒質よりも明るく見える。 同様に、負相関の場合は粒子同士の重なり合いは小さくなり、 減衰は指数的なものよりも大きくなり、ボリュームはより暗く見える。

これをモデル化する上で厄介なのが、 ある粒子の位置が別の粒子の位置に対して情報量を持つことである。 このとき、 光が媒質内で減衰する確率はこの光の過去の状態 (経路) に依存し、記憶を持つと言える (memory effect)。 今回は特に、光がボリューム内で最後に反射した位置から進んだ距離を t とすると、 減衰確率は t の関数として表さなければならない。

これをモデル化するために消散確率を拡張したものが、消散確率微分 (differential probability of extinction) であり、 これは下の式で表される。

 \displaystyle \Sigma(t) = \frac{ p(t) }{ 1 - \int_0^t p(s) ds }

 (p(t) = \mu exp(-\mu t))

 p(t) は距離 t で消散が起こる確率密度であり、 要するにこれは『t までに消散されなかった場合に、t で消散される確率密度』を表す条件付き確率密度である。

もちろんこれは粒子の空間的相関性が無い場合には通常の消散係数  \mu に一致する。

The Generalized Boltzmann Equation (GBE)

粒子の空間的相関性を扱えるようにするため、上記の消散確率微分  \Sigma(t) を RTE に導入する。  \Sigma(t) は t の関数であるため、このとき輝度も t の関数  L(\boldsymbol{x}, \omega_o, t) として拡張する。  L(\boldsymbol{x}, \omega_o, t) \boldsymbol{x} \omega_o 方向で到達する直前に 最後の拡散イベントからちょうど t だけ直進していた、という光の輝度を表す。 単位は輝度の t に関する微分となる。

RTE に適用すると以下のようになる。

 \displaystyle \frac{d}{dt} L(\boldsymbol{x}, \omega_o, t) + \omega_o \cdot \nabla L(\boldsymbol{x}, \omega_o, t) + \Sigma(t) L(\boldsymbol{x}, \omega_o, t) = 0

 \displaystyle L(\boldsymbol{x},\omega_o,0) = \int_0^{\infty} \Sigma_s(t) \int_\Omega L(\boldsymbol{x},\omega_o,t) f_r(\omega_i,\omega_o) d\omega_i dt + Q(\boldsymbol{x},\omega_o)

( \Sigma_s は消散確率微分  \Sigma(t)アルベド \Lambda の積)

二個目の式は、点 x における入射光の拡散と発光についてを表している。

一個目の式は、 第一項が輝度の微分微分になっているのでややこしいが、 要するに

(拡散イベント位置についての L の微分) + (x についての L の微分) + (x においての L の減衰) = 0

つまり、

-(入ってくる L) + (出ていく L) + (L の減衰) = 0

となっており、当然だが形自体は RTE と似通っている。 (t の負方向が  \omega_o を向くので第一項の符号は上記のようになるのに注意)。

この拡張された RTE を、The Generalized Boltzmann Equation (GBE) と呼ぶ。  \Sigma(t) を定数とすれば、当然 GBE は RTE と同じになり、 また  \Sigma(t) を方向の関数に拡張すれば異方性媒質の表現が可能になる。

GBE のリミテーション

だが GBE はいくつか単純化した仮定を置いており、 そのままボリュームレンダリングに適用することはできない。 GBE が前提とする仮定は以下である。

  • 媒質は均一で、無限に広がっており境界が無い。
  • 位相関数とアルベドは t に依存していない。つまり、異なる空間的相関性を持つ複数種類の粒子が混在している場合を扱えない。(粒子の空間的相関性が異なる場合、t が小さいうちは粒子 A の影響が強く、t が大きくなると粒子 B の影響も大きくなる、などの現象が起こる)
  • ボリュームの発する光も光源となる粒子から発されていると考えられるが、このとき、光源の粒子は光を散乱する通常の粒子と同じ空間的相関性を持っていると考える。これは多くのケースでは正しくない。

提案手法

拡張 GBE

上記の 1,2 つ目の GBE のリミテーションを解決するため、 消散確率微分  \Sigma に加えて位相関数  f_rアルベド  \Lambda を t の関数とする (以下、簡単のため  B(\boldsymbol{x},\omega_i,\omega_ot) = \Lambda(\boldsymbol{x},t) \Sigma(\boldsymbol{x},t) f_r(\boldsymbol{x},\omega_i,\omega_o,t) とする)。 また、3 つ目のリミテーションに対処するため、 ボリュームからの発光は離散的な光源の和で表現できると仮定し、 L を他から x に入射して拡散される光と x から発される光とに分けて

 \displaystyle  L(\boldsymbol{x}, \omega_o, t) = L_S(\boldsymbol{x}, \omega_o, t) + \sum_j L_{Q_j}(\boldsymbol{x}, \omega_o, t)

という表現ができるとする( L_S(\boldsymbol{x}, \omega_o, t) は拡散された輝度、 L_{Q_j}(\boldsymbol{x}, \omega_o, t) はボリューム内の光源  Q_j から発された輝度)。

以上を GBE に適用すると以下の式が得られる。

 \displaystyle \frac{d}{dt} L(\boldsymbol{x}, \omega_o, t) + \omega_o \cdot \nabla L(\boldsymbol{x}, \omega_o, t) + \Sigma_S(\boldsymbol{x},t) L_S(\boldsymbol{x}, \omega_o, t) + \sum_j \Sigma_{Q_j}(\boldsymbol{x},t) L_{Q_j}(\boldsymbol{x}, \omega_o, t) = 0

 \displaystyle L_S(\boldsymbol{x},\omega_o,0) = \int_0^{\infty} \int_\Omega \Bigl( B_S(\boldsymbol{x},\omega_i,\omega_o,t) L_S(\boldsymbol{x},\omega_o,t) + \sum_j
B_{Q_j}(\boldsymbol{x},\omega_i,\omega_o,t) L_{Q_j}(\boldsymbol{x},\omega_i,t) \Bigr) d\omega_i dt

 L_{Q_j} (\boldsymbol{x},\omega_o,0) = Q_j(\boldsymbol{x},\omega_o)

( \Sigma_S,  \Sigma_{Q_j} は拡散された光、また発された光それぞれの消散確率微分。)

上式を用いると、 媒質の境界面での光のやりとりは 媒質境界面に届いたレイを光源  Q_j として扱ってやったりといった 境界条件を定めることによって可能になる。

混合媒質について

アルベドと位相関数は t の関数として拡張された。 これは、複数種の粒子が混在している媒質を扱うための拡張なので、 消散確率微分も含めて混合媒質についての定式化を以下で行う。 (粒子が 1 種類の場合は、アルベドと位相関数は t についての関数にする必要はない。)

まず、 媒質を構成する粒子種類の集合を  \mathcal{P} とすると、 媒質の消散確率微分  \Sigma(\boldsymbol{x},t) は個々の粒子種類  k \in \mathcal{P} の 消散確率微分  \Sigma_k(\boldsymbol{x},t) の重み付き和として、

 \displaystyle \Sigma(\boldsymbol{x},t) = \sum_{k \in \mathcal{P}} w_k \Sigma_k(\boldsymbol{x},t)

のように定義できる。

このときアルベドは、同様に各粒子種類のアルベドの重み付き和となるが、 重みとしてさらに消散確率微分を考慮して

 \displaystyle \Lambda(\boldsymbol{x},t) = \sum_{k \in \mathcal{P}} \frac{w_k \Sigma_k(\boldsymbol{x},t)}{\sum_{k \in \mathcal{P}} w_k \Sigma_k(\boldsymbol{x},t)} \Lambda_k(\boldsymbol{x},t)

となる。

位相関数は、重みとしてさらにアルベドを考慮して

 \displaystyle f_r(\boldsymbol{x}, \omega_i, \omega_o, t) = \sum_{k \in \mathcal{P}} \frac{w_k \Sigma_k(\boldsymbol{x},t) \Lambda_k(\boldsymbol{x},t)}{\sum_{k \in \mathcal{P}} w_k \Sigma_k(\boldsymbol{x},t)\Lambda_k(\boldsymbol{x},t)} f_{r,k}(\boldsymbol{x}, \omega_i, \omega_o, t)

となる。

消散確率微分関数の定式化

残った問題は、 消散確率微分  \Sigma(t) をどのようにモデル化するかということである。 以下では  \Sigma(t) を透過率  T(t) から求める。

まず、 ある点に入射してくる輝度  L_i が確率密度  p_L(L_i) に従って現れる確率的な変数だと考える。 このとき、消散係数  \mu を考えることができ、これは個々のレイが通った経路に応じて決まる確率的な変数になる。その確率密度関数 p_\tau(\mu) とする。

話を簡単化するため 媒質が均一で粒子間に空間的正相関を持っており、 さらに  p_L(L_i) p_\tau(\mu) の間に相関性が無いとする。

すると、途中の式は省くが、 t だけ直進した光の透過率は

 \displaystyle T(t) = \int_0^\infty e^{-\mu t} p_\tau(\mu) d\mu

として表すことができる。 これは、空間的相関性のために個々の光路に応じて  \mu が変化するが、 その各  \mu についての減衰率  e^{-\mu t} の期待値が透過率になる、と解釈できる。多分。

ここで、 p(t) = \left| \frac{dT(t)}{dt} \right| \Sigma(t) = p(t) / T(t) という関係から、

 \displaystyle \Sigma(t) = \frac{\int_0^\infty \mu e^{-\mu t} p_\tau(\mu) d\mu}{\int_0^\infty e^{-\mu t} p_\tau(\mu) d\mu}

が得られる。 (なお、  p_L(L_i) p_\tau(\mu) の間に相関性が無いと仮定しているため、 光源  Q_j と粒子の間には相関性がなく、 結局  \Sigma(t) = \Sigma_S(t) = \Sigma_{Q_j}(t) となる。)

以上から  \Sigma(t) のモデル化は  p_\tau(\mu) のモデル化に帰着されたことになった。 ここでさらに、 粒子の平均断面積を  \sigma とすると、 消散係数  \mu (体積あたりの粒子断面積総和) と 粒子の濃度  \mathcal{C} (体積あたりの粒子個数) の間の関係  \mu p_\tau(\mu) = \mathcal{C} p_{\mathcal{C}}(\mathcal{C}) \sigma から、 濃度の確率密度  p_\mathcal{C}(\mathcal{C}) のモデル化へと帰着でき、 この  p_\mathcal{C}(\mathcal{C}) をガンマ分布として定義するといい感じになるというのが実験的に確かめられているらしい。 またそうすることで、 \Sigma(t) も解析的に解けることになる。

結果とか

ページ頭にも載せたが、 無相関、正相関、負相関についてのレンダリングは以下のようになる。 任意の負相関についての消散確率微分のモデル化については解決されていないが、 完全な負相関の場合は線形に透過率が線形に下がることにより一応解けるためそれを使っているらしい。 f:id:mizuooon:20180605162830p:plain

無相関媒質のパラメータをいじって正相関の見た目に頑張って近づけようとしてみたものが以下。 透過のニュアンスが無相関媒質では再現しきれていないのがわかる。 f:id:mizuooon:20180612133605p:plain

相関についての計算を行う追加コストはほとんど無視できるらしい。

今後の課題としては、 負相関の消散確率微分のモデル化であったり、 不均質な連続媒体のサポートなどが挙げられている。