2017年8月13日日曜日

Layered Diffuse Surface (Interfaced Lambertian BRDF)

Layered Diffuse Surface (Interfaced Lambertian BRDF)

このページは レイトレ合宿5‽ アドベントカレンダー第9週目の記事です.

本ページでは,モンテカルロレイトレーシング法において, 上図のような鏡面反射する層を持った拡散反射基板の表面 (Interfaced Lambertian BRDF) をレンダリングする方法を説明します. Interfaced Lambertian BRDF の実装に必要な パラメータや数式,重点的サンプリングの方法を説明します.

記号について

本ページでは以下の記号は定義済として扱います.

記号 説明
\(\omega_{n}\) Macrosurfaceの単位法線ベクトル
\(\omega_{m}\) Microsurfaceの単位法線ベクトル,本ページではマイクロファセット法線と呼ぶ
\(\omega_{i}\) 入射方向の単位ベクトル
\(\omega_{o}\) 反射方向の単位ベクトル
\(f_{r}\) BRDF
\(n_{i}\) 入射側の媒体の屈折率
\(n_{o}\) 透過側の媒体の屈折率
\(n\) \(\frac{n_{o}}{n_{i}}\)
\(F\) フレネル項

本モデルについて

Interfaced Lambertian BRDF は,下図のようにフレネル反射層と拡散反射基板から構成されます. このモデルは,フレネル反射層の境界での鏡面反射と 層内部の多重反射後の散乱を考慮しています. フレネル反射層と拡散反射基板はそれぞれ鏡面反射,拡散反射のMicrofacetモデルで表現され, 必要なパラメータもそれぞれのMicrofacetモデルのパラメータを組み合わせるだけの シンプルなモデルとなっています.

Interfaced Lambertian BRDF

Interfaced Lambertian BRDF の評価は, 入射光に対する鏡面反射 \(f_{s}\) と, 多重反射後の散乱 \(f_{b}\) に分けて行います.

\[ \begin{align} f_{r} \left( x, \omega_{i}, \omega_{o} \right) = f_{s} \left( x, \omega_{i}, \omega_{o} \right) + f_{b} \left( x, \omega_{i}, \omega_{o} \right) \tag{1} \end{align} \]

入射光に対する鏡面反射 \(f_{s}\) には, 鏡面のMicrofacetモデルを用います (Beckman,GGXなど) .

\[ \begin{align} f_{s} \left( x, \omega_{i}, \omega_{o} \right) = \frac{1}{4 \left( \omega_{i} \cdot \omega_{n} \right) \left( \omega_{o} \cdot \omega_{n} \right)} F \left( \omega_{i}, \omega_{m} \right) G_{2} \left( \omega_{i}, \omega_{o}, \omega_{m} \right) D \left( \omega_{m} \right) \tag{2} \end{align} \]

多重反射後の散乱 \(f_{b}\) は, 拡散反射面のMicrofacetモデルに, 入射光の内部への屈折 \(T \left( \omega_{i}, \omega_{m} \right)\) と 内部での多重反射 \(\frac{1}{\pi n^{2} \left( 1 - r_{i} \rho \right) }\) と, 内部から外部へ散乱する時の屈折 \(T \left( \omega_{o}, \omega_{m} \right)\) を 考慮したものになります.

\[ \begin{align} f_{b} \left( x, \omega_{i}, \omega_{o} \right) = & \frac{\rho}{\pi n^{2} \left( 1 - r_{i} \rho \right)} \int T \left( \omega_{i}, \omega_{m_{b}} \right) T \left( \omega_{i}, \omega_{m_{b}} \right) D \left( \omega_{m_{b}} \right) G \left( \omega_{i}, \omega_{o}, \omega_{m_{b}} \right) \frac{ \left| \omega_{i} \cdot \omega_{m_{b}} \right| \left| \omega_{o} \cdot \omega_{m_{b}} \right| } { \left| \omega_{i} \cdot \omega_{n} \right| \left| \omega_{o} \cdot \omega_{n} \right| } d \omega_{m_{b}} \\ r_{i} = & 1 - \frac{1 - r_{e}}{n^{2}} \\ r_{e} = & \frac{1}{2} - \frac{2 n^{3} \left( n^{2} + 2 n - 1 \right)} {\left( n^{2} + 1 \right) \left( n^{4} - 1 \right)} + \frac{\left( n - 1 \right) \left( 3 n + 1 \right)} {6 \left( n + 1 \right)^{2}} + \frac{8 n^{4} \left( n^{4} + 1 \right)} {\left( n^{2} + 1 \right) \left( n^{4} - 1 \right)^{2}} ln \left( n \right) + \frac{n^{2} \left( n^{2} - 1 \right)^{2}} {\left( n^{2} + 1 \right)^{3}} ln \left( \frac{n - 1}{n + 1} \right) \tag{3} \end{align} \]

この式には積分が含まれているため,モンテカルロ法を用いて Microfacet法線 \(\omega_{m_{b}}\) をサンプリングして解きます.

\[ \begin{align} f_{b} \left( x, \omega_{i}, \omega_{o} \right) = \frac{\rho}{\pi n^{2} \left( 1 - r_{i} \rho \right)} \frac{ T \left( \omega_{i}, \omega_{m_{b}} \right) T \left( \omega_{i}, \omega_{m_{b}} \right) D \left( \omega_{m_{b}} \right) G \left( \omega_{i}, \omega_{o}, \omega_{m_{b}} \right) \frac{ \left| \omega_{i} \cdot \omega_{m_{b}} \right| \left| \omega_{o} \cdot \omega_{m_{b}} \right| } { \left| \omega_{i} \cdot \omega_{n} \right| \left| \omega_{o} \cdot \omega_{n} \right| }} {p \left( \omega_{m_{b}} \right) } \tag{4} \end{align} \]

重点的サンプリング (Importance Sampling)

Interfaced Lambertian BRDF の重点的サンプリングは, フレネル反射層の反射分布,もしくは拡散反射基板の反射分布に従ってサンプリングを行います.

まず最初は,フレネル反射層,拡散反射基板どちらの反射分布に従ってサンプリングを行うかを決定します. ここでは,フレネル反射層,拡散反射基板それぞれの総反射率をもとに確率的に決定します.

\[ \begin{align} R_{s} = & r_{e} \\ R_{b} = & \frac{ \left( 1 - r_{e} \right) ^{2} \rho } { n^{2} \left( 1 - r_{i} \rho \right) } \tag{5} \end{align} \]

\(R_{s}\)\(R_{b}\) はそれぞれフレネル反射層,拡散反射基板の総反射率です. \(u_{1} \in [0, 1)\) となる一様乱数を用いて, \(u_{1} < \frac{R_{s}}{ \left( R_{s} + R_{b} \right) }\) となる場合は フレネル反射層の反射分布に従ったサンプリングを行い,そうでない場合は, 拡散反射基板の反射分布に従ったサンプリングを行います.

フレネル反射層の反射分布のpdf \(p_{s} \left( \omega_{o} \right)\) や 反射方向のサンプリング方法は層のMicrofacetモデルに従う (Beckman,GGXなど). 反射方向のサンプリングの方法については,アドベントカレンダー第8回Rough Specular Surface (GGX BRDF) を参照してください.

拡散反射基板の反射分布のpdf \(p_{b} \left( \omega_{o} \right)\) は,

\[ \begin{align} p_{b} \left( \omega_{o} \right) = \frac{ \left( \omega_{o} \cdot \omega_{n} \right) }{\pi} \tag{6} \end{align} \]

であり, \(p_{b} \left( \omega_{o} \right)\) に従って反射方向をサンプリングします. 反射方向のサンプリング方法は Smooth Diffuse Surface (Lambert BRDF) を参照してください.

最終的な反射方向のpdf \(p \left( \omega_{o} \right)\) は,

\[ \begin{align} p \left( \omega_{o} \right) = \frac{R_{s}}{ \left( R_{s} + R_{b} \right) } p_{s} \left( \omega_{o} \right) + \frac{R_{b}}{ \left( R_{s} + R_{b} \right) } p_{b} \left( \omega_{o} \right) \tag{7} \end{align} \]

となります.

まとめ

パラメータ

名称 記号
拡散反射基板の反射率 \(\rho\) \(0 \le \rho \le 1\)
表面粗さ \(roughness\) \(0 < roughness \le 1\)

レンダリング

名称 記号
BRDF \(f_{r} \left( x, \omega_{i}, \omega_{o} \right)\) \(f_{s} \left( x, \omega_{i}, \omega_{o} \right) + f_{b} \left( x, \omega_{i}, \omega_{o} \right)\)
pdf \(p \left( \omega_{o} \right)\) \(\frac{R_{s}}{ \left( R_{s} + R_{b} \right) } p_{s} \left( \omega_{o} \right) + \frac{R_{b}}{ \left( R_{s} + R_{b} \right) } p_{b} \left( \omega_{o} \right)\)
反射レイのWeight \(\frac{f_{r} \left( x, \omega_{i}, \omega_{o} \right) \left( \omega_{o} \cdot \omega_{n} \right)} {p \left( \omega_{o} \right)}\) \(\frac{f_{r} \left( x, \omega_{i}, \omega_{o} \right) \left( \omega_{o} \cdot \omega_{n} \right)} {p \left( \omega_{o} \right)}\) (特にキャンセルできる項は無い)

以上,Interfaced Lambertian BRDFについて簡単にまとめました.

参考文献

  1. Daniel Meneveaux: Rendering Rough Opaque Materials with Interfaced Lambertian Microfacets (2017)

2016年8月7日日曜日

BRDFのテスト

BRDFのテスト

このページは レイトレ合宿4!? アドベントカレンダー第7週目の記事です. レイトレ合宿まで1ヶ月をきりました. そろそろ提出用のシーン作成やレンダラーのデバッグに集中しているころではないでしょうか.

私はレンダラーの実装に不具合があるかどうかは, レンダリングした画像を見て違和感を持つかどうかで判断しています. この方法は主観的であり, 画像に違和感を持たなくても不具合が含まれている可能性は大いにあります. そこで,自動テストを実装して客観的に不具合の判定をするようにします.

このページでは物理ベースなBRDFの自動テストについて考えます.

物理ベースなBRDF

初めに,物理ベースなBRDFについて確認しておきます.

定義済み記号 説明
\(\omega_{n}\) 法線ベクトル
\(\omega_{i}\) 入射方向の単位ベクトル
\(\omega_{o}\) 反射方向の単位ベクトル
\(f_{r}\) BRDF

物理ベースなBRDFは,反射する放射エネルギーが負になることはないので

\[ \begin{align} \forall \boldsymbol{\omega_{i}} \ \forall \boldsymbol{\omega_{o}}, \ f_{r}(x, \boldsymbol{\omega_{o}} \leftarrow \boldsymbol{\omega_{i}}) \ cos \theta_{\boldsymbol{\omega_{o}}} \ge 0 \end{align} \]

\(cos \theta_{\boldsymbol{\omega_{o}}} \ge 0\) であるため,

\[ \begin{align} \forall \boldsymbol{\omega_{i}} \ \forall \boldsymbol{\omega_{o}}, \ f_{r}(x, \boldsymbol{\omega_{o}} \leftarrow \boldsymbol{\omega_{i}}) \ge 0 \tag{1} \end{align} \]

となります. また,物理ベースなBRDFは以下の2つの特性を持ちます.

ヘルムホルツの相反性 (Helmholtz Reciprocity)

BRDFの値は入射・出射方向が入れ替わっても変化しません.

\[ \begin{align} f_{r}(x, \boldsymbol{\omega_{o}} \leftarrow \boldsymbol{\omega_{i}}) = f_{r}(x, \boldsymbol{\omega_{i}} \leftarrow \boldsymbol{\omega_{o}}) \tag{2} \end{align} \]

エネルギー保存則 (Energy Conservation)

ある点から反射した総出射光束は総入射光束以下となります.

\[ \begin{align} \forall \boldsymbol{\omega_{i}}, \ \int_{\Omega} f_{r}(x, \boldsymbol{\omega_{o}} \leftarrow \boldsymbol{\omega_{i}}) \ cos \theta_{\boldsymbol{\omega_{o}}} d \sigma(\boldsymbol{\omega_{o}}) \le 1 \tag{3} \end{align} \]

モンテカルロ法によるテスト

\((1)\)\((2)\)\((3)\)の自動テストを実装します. 式\((3)\)は解析的に解くことが難しいため,モンテカルロ法を用いてテストを行います. モンテカルロ法を用いた場合,式\((3)\)は以下のように書きかえられます.

\[ \begin{align} \forall \boldsymbol{\omega_{i}}, \frac{1}{N} \sum_{k=1}^{N} \frac{f_{r}(x, \boldsymbol{\omega_{ok}} \leftarrow \boldsymbol{\omega_{i}}) cos \theta_{\boldsymbol{\omega_{ok}}}} {p_{u} \left( \boldsymbol{\omega_{ok}} \right)} \le 1 \tag{4} \end{align} \]

ここで,\({p_{u} \left( \boldsymbol{\omega_{o}} \right)}\)\(\boldsymbol{\omega_{o}}\)の確率密度関数 (PDF) です. このテストでは,半球上の方向を一様にサンプルするようにします \(\left({p_{u} \left( \boldsymbol{\omega_{o}} \right)} = \frac{1}{2 \pi} \right)\)

C++言語で GoogleTest を用いたテストの擬似コードを記述します.

// エネルギー保存則のテスト
TEST(BrdfTest, EnergyConservationTest)
{
  Brdf brdf; // テスト対象であるBRDFの宣言,初期化
  const int N = 1000000;
  // 様々な入射方向に対してテストを行う
  for (int i = 0; i < N; ++i) {
    const float energySum = 0.0;
    const Vector3 vin = sampleDirectionUniformly(); // 半球上の方向を一様にサンプル
    for (int o = 0; o < N; ++i) {
      const Vector3 vout = sampleDirectionUniformly();
      const float f = brdf.f(vin, vout); // BRDFの評価
      ASSERT_LE(0.0, f); // 式(1)のテスト
      const float cosTheta = dot(normal, vout);
      energySum += f * cosTheta * (2.0 * pi); // pi は円周率
    }
    const float energy = energySum / (float)N;
    ASSERT_GE(1.0, energy); // 式(4)のテスト
  }
}
// ヘルムホルツの相反性のテスト
TEST(BrdfTest, HelmHoltzReciprocityTest)
{
  Brdf brdf;
  const int N = 1000000;
  for (int i = 0; i < N; ++i) {
    const Vector3 vin = sampleDirectionUniformly();
    const Vector3 vout = sampleDirectionUniformly();
    const float f1 = brdf.f(vin, vout); // BRDFの評価
    const float f2 = brdf.f(-vout, -vin); // 方向を入れ替えたBRDFの評価
    const float error = 1.0e-5; // 許容できる誤差を定義
    ASSERT_FLOAT_EQ(f1, f2, error); // 式(2)のテスト
  }
}

実際にテストを実装する場合は, モンテカルロ法による\(N\)回試行の和を求める際に情報落ちの危険があるため, 補正加算も考慮したほうがいいと思います.

重点的サンプリング (Importance sampling) のテスト

BRDFによっては,反射方向を重点的サンプリングする式が定義してあることがあります. 重点的サンプリングの際の反射方向のPDFを \(p \left( \boldsymbol{\omega_{o}} \right)\) とします. 確率密度関数の性質から,

\[ \begin{align} \int_{\Omega} p \left( \boldsymbol{\omega_{o}} \right) d \sigma \left( \boldsymbol{\omega_{o}} \right) = 1 \tag{5} \end{align} \]

となります.式\((3)\)と同様に,解析的に解くことが難しいためモンテカルロ法を 用いて以下のように書きかえます.

\[ \begin{align} \frac{1}{N} \sum_{k=1}^{N} \frac{p \left( \boldsymbol{\omega_{ok}} \right)} {p_{u} \left( \boldsymbol{\omega_{ok}} \right)} \approx 1 \tag{6} \end{align} \]

擬似コードを以下に記述します.

// PDFのテスト
TEST(BrdfTest, PdfTest)
{
  Brdf brdf;
  const int N = 1000000;
  // 様々な入射方向に対してテストを行う
  for (int i = 0; i < N; ++i) {
    float pdfSum = 0.0;
    const Vector3 vin = sampleDirectionUniformly();
    for (int o = 0; o < N; ++o) {
      const Vector3 vout = sampleDirectionUniformly();
      const float pdf = brdf.pdf(vin, vout); // 反射方向のPDFを評価
      pdfSum += pdf * (2.0 * pi);
    }
    const float pdf = pdfSum / (float)N;
    const float error = 1.0e-5; // 許容できる誤差を定義
    ASSERT_FLOAT_NEAR(1.0, pdf, error); // 式(6)のテスト
  }
}

終わりに

BRDFのテストについて考えました. このテストで確認できたことは,物理に基づいたBRDFの定義を満たしているか ということであり,BRDFの反射分布が正しいかというところまでは確認できていません. ここからテスト項目を増やしていって,確認できることを広げていこうと思います.

2016年3月5日土曜日

Rough Specular Surface (GGX BRDF)

Rough Specular Surface (GGX BRDF)

本ページでは,モンテカルロレイトレーシング法において, 金属など絶縁体の粗い表面を表現する方法を説明する. 言い換えると,粗い表面モデルのBRDFの実装に必要なパラメータや数式, 重点的サンプリングの方法を説明する. 本ページでは,粗い表面モデルとして[1]のGGXモデルを用いる.

記号について

本ページでは以下の記号は定義済として扱う.

記号 説明
\(\omega_{n}\) Macrosurfaceの単位法線ベクトル
\(\omega_{m}\) Microsurfaceの単位法線ベクトル,本ページではマイクロファセット法線と呼ぶ
\(\omega_{i}\) 入射方向の単位ベクトル
\(\omega_{o}\) 反射方向の単位ベクトル
\(f_{r}\) BRDF
\(n_{i}\) 入射側の媒体の屈折率
\(n_{o}\) 透過側の媒体の屈折率
\(\eta_{o}\) 透過側の媒体の消衰係数
\(n\) \(\frac{n_{o}}{n_{i}}\)
\(\eta\) \(\frac{\eta_{o}}{n_{i}}\)
\(F\) フレネル項
\(\chi^{+} \left( x \right)\) \(\begin{cases} 1 & (0 < x) \\ 0 & (otherwise) \end{cases}\)

マイクロファセットモデルとGGX

マイクロファセットモデルは様々な 方向を向いた微小平面で構成される表面(Microsurface)を 単純な平面(Macrosurface)で近似したモデルである. このモデルでは,単純化のため波動光学や表面での2回以上の反射は無視する. マイクロファセットモデルは2つの統計量, 法線分布関数 \(D\) と Shadowing-Masking Function \(G\) を用いて表される.

表面の粗さについて

GGXは微小平面で見れば鏡面反射しかしていないが, 様々な法線を持った微小平面から構成されているため 表面が粗いほど拡散反射をしているように見える. 表面の粗さは \(roughness\) パラメータを用いて制御する. \(roughness\) がとる値の範囲は \(0 < roughness \le 1\) である. 値が \(0\) に近いほど理想的な鏡面反射に近くなり, \(1\) に近いほど理想的な拡散反射に近くなる.

また,計算時に用いる表面粗さの値は, \(roughness\) を直接用いるのではなく, [3]で解説してあるように

\[ \begin{align} \alpha = \left( roughness \right)^{2} \tag{1} \end{align} \]

を用いるほうがより視覚的に線形比例した変化となる. 本ページでも, \(\alpha\) を計算時の表面粗さの値として用いる.

法線分布関数 (Normal Distribution Function) \(D\)

法線分布関数 \(D\) はMicrosurface上にできる法線の統計分布を表したものである. GGXでは \(D\) は以下の式で表される.

\[ \begin{align} D \left( \omega_{m} \right) = \chi^{+} \left( \omega_{n} \cdot \omega_{m} \right) \frac{\alpha^{2}} {\pi \left( \left( \alpha^{2} - 1 \right) \left( \omega_{n} \cdot \omega_{m} \right)^{2} + 1 \right)^{2}} \tag{2} \end{align} \]

ここでは簡易表記のため, \(D \left( \omega_{n}, \omega_{m}, \alpha \right)\)\(D \left( \omega_{m} \right)\) と表記している.

Shadowing-Masking Function \(G\)

Shadowing-Masking Function は[2]の論文では V-cavityモデルとSmithモデルの2通りのモデルを載せているが, ここではSmithモデルを説明する(V-cavityモデルについては[2]を参照).

Shadowing-Masking Function \(G\) はある法線を持つMicrosurfaceが可視かどうかの比率 を表したものである. \(G_{2}\)は入射方向,反射方向の双方向から可視かどうかの比率を表しており, \(G_{1}\)は入射方向もしくは反射方向のどちらか単方向から可視かどうかの比率を表している.

\[ \begin{align} G_{2} \left( \omega_{i}, \omega_{o}, \omega_{m} \right) = & G_{1} \left( \omega_{i}, \omega_{m} \right) G_{1} \left( \omega_{o}, \omega_{m} \right) \tag{3} \\ G_{1} \left( \omega, \omega_{m} \right) = & \chi^{+} \left( \left( \omega \cdot \omega_{n} \right) \left( \omega \cdot \omega_{m} \right) \right) \frac{2 \left( \omega \cdot \omega_{n} \right)} { \left( \omega \cdot \omega_{n} \right) \sqrt{\alpha^{2} + \left( 1 - \alpha^{2} \right) \left( \omega \cdot \omega_{n} \right)^{2} }} \tag{4} \end{align} \]

ここでは簡易表記のため, \(G_{1} \left( \omega, \omega_{n}, \omega_{m}, \alpha \right)\)\(G_{1} \left( \omega, \omega_{m} \right)\) と表記している.

GGX BRDF

GGXのBRDFは以下の式になる.

\[ \begin{align} f_{r} \left( x, \omega_{i}, \omega_{o} \right) = \frac{1}{4 \left( \omega_{i} \cdot \omega_{n} \right) \left( \omega_{o} \cdot \omega_{n} \right)} F \left( \omega_{i}, \omega_{m} \right) G_{2} \left( \omega_{i}, \omega_{o}, \omega_{m} \right) D \left( \omega_{m} \right) \tag{5} \end{align} \]

ここでは簡易表現のため, \(F \left( \omega, \omega_{m}, n, \eta \right)\)\(F \left( \omega, \omega_{m} \right)\) と表記している.

ExplicitConnection のように, マイクロファセット法線 \(\omega_{m}\) がわからない状態で \(f_{r} \left( x, \omega_{i}, \omega_{o} \right)\) を評価する場合は, \(\omega_{m}\) は 反射ハーフベクトル \(\omega_{h_{r}}\) を用いて,

\[ \begin{align} \omega_{m} = & \frac{\omega_{h_{r}}} {\left| \left| \omega_{h_{r}} \right| \right|} \tag{6} \\ \omega_{h_{r}} = & \omega_{i} + \omega_{o} \tag{7} \end{align} \]

から求める.

重点的サンプリング (Importance Sampling)

反射方向 \(\omega_{o}\) の重点的サンプリングを行う. GGX BRDF では, マイクロファセット法線 \(\omega_{m}\) を重点的サンプリングし, 入射方向 \(\omega_{i}\)\(\omega_{m}\) から正反射方向 \(\omega_{o}\) を求める.

\(\omega_{m}\) は 投影法線分布関数(Projected Normal Distribution Function) \(D_{\omega_{i}} \left( \omega_{m} \right)\) に従ってサンプリングする. \(D_{\omega_{i}} \left( \omega_{m} \right)\) は, \(D \left( \omega_{m} \right)\) の中で 入射方向 \(\omega_{i}\) から見える法線のみの分布を示したもので,

\[ \begin{align} D_{\omega_{i}} \left( \omega_{m} \right) = \frac{ \left( \omega_{i} \cdot \omega_{m} \right) } { \left( \omega_{i} \cdot \omega_{n} \right) } G_{1} \left( \omega_{i}, \omega_{m} \right) D \left( \omega_{m} \right) \tag{8} \end{align} \]

で定義される (ここでは簡易表記のため, \(D_{\omega_{i}} \left( \omega_{n}, \omega_{m} \right)\)\(D_{\omega_{i}} \left( \omega_{m} \right)\) と表記している). マイクロファセット法線 \(\omega_{m}\) の pdf \(p \left( \omega_{m} \right) = D_{\omega_{i}} \left(\omega_{m} \right)\) から, 反射方向 \(\omega_{o}\) とそのpdf \(p \left( \omega_{o} \right)\)

\[ \begin{align} \omega_{o} = & 2 \left( \omega_{i} \cdot \omega_{m} \right) \omega_{m} - \omega_{i} \tag{9} \\ p \left( \omega_{o} \right) = & \left| \left| \frac{\delta \omega_{m}}{\delta \omega_{i}} \right| \right| p \left( \omega_{m} \right) \\ = & \frac{1}{4 \left( \omega_{i} \cdot \omega_{m} \right)} D_{\omega_{i}} \left( \omega_{m} \right) \tag{10} \end{align} \]

で求められる. ここで,\(\left| \left| \frac{\delta \omega_{m}}{\delta \omega_{i}} \right| \right| = \frac{1}{4 \left( \omega_{m} \cdot \omega_{i} \right)}\) はヤコビアンである.

重点的サンプリングの際の反射レイのWeightは

\[ \begin{align} \frac{f_{r} \left( x, \omega_{i}, \omega_{o} \right) \left( \omega_{o} \cdot \omega_{n} \right)} {p \left( \omega_{o} \right)} = F \left( \omega_{i}, \omega_{m} \right) G_{1} \left( \omega_{o}, \omega_{m} \right) \tag{11} \end{align} \]

となる.

\(\omega_{m}\) のサンプリング

\(D_{\omega_{i}}(\omega_{m})\) に従ってマイクロファセット法線 \(\omega_{m}\) をサンプリングする. サンプリング方法はV-cavity法とSmith法の2通りの方法があるが, ここではSmith法を載せる(V-cavity法については[2]を参照).

注意として,このサンプリング方法は macrosurfaceの法線が \(\omega_{n} = (0, 0, 1)\) の時の 入射角 \(\omega_{i} = (x_{i}, y_{i}, z_{i})\) を用いる. \(\omega_{n}\) がそれ以外の場合は \(\omega_{n} = (0, 0, 1)\) となるように \(\omega_{i}\) の基底変換を行う. \(\omega_{m}\) のサンプリング後は \(\omega_{i}\)\(\omega_{m}\) を 元の基底に変換する.

Smith法による\(\omega_{m}\)のサンプリング方法は 大きく分けて以下の5つのプロセスからなる.

  1. Stretch \(\omega_{i}\)
  2. Sample \(P_{\omega_{i}}^{22}(x_{\tilde{m}}, y_{\tilde{m}}, 1)\)
  3. Rotate \(\omega_{m}\)
  4. Unstretch \(\omega_{m}\)
  5. Normalize \(\omega_{m}\)

Stretch \(\omega_{i}\)

\[ \begin{align} \omega'_{i} \equiv \frac{(\alpha x_{i}, \alpha y_{i}, z_{i})} {\sqrt{(\alpha x_{i})^{2} + (\alpha y_{i})^{2} + z_{i}^2}} \end{align} \]

Sample \(P_{\omega_{i}}^{22}(x_{\tilde{m}}, y_{\tilde{m}}, 1)\)

\(x_{\tilde{m}}\) のサンプリング

\[ \begin{align} A \equiv & u_{1} \left( 1 + \sqrt{1 + tan^{2} \theta_{\omega'_{i}}} \right) - 1 \\ B \equiv & tan \theta_{\omega'_{i}} \\ x_{\tilde{m}1} \equiv & \frac{B}{A^{2} - 1} - \sqrt{\left(\frac{B}{A^{2} - 1}\right)^{2} - \frac{A^{2} - B^{2}}{A^{2} - 1}} \\ x_{\tilde{m}2} \equiv & \frac{B}{A^{2} - 1} + \sqrt{\left(\frac{B}{A^{2} - 1}\right)^{2} - \frac{A^{2} - B^{2}}{A^{2} - 1}} \\ x_{\tilde{m}} \equiv & \begin{cases} x_{\tilde{m}1} & (A < 0 \ \ or \ \ cot \theta_{\omega'_{i}} < x_{\tilde{m}2}) \\ x_{\tilde{m}2} & (otherwise) \end{cases} \end{align} \]

ここで,\(u_{1}\)\([0, 1)\) の乱数である.

\(y_{\tilde{m}}\) のサンプリング

\[ \begin{align} s \equiv & \begin{cases} 1 & (0.5 \lt u_{2}) \\ -1 & (otherwise) \end{cases} \\ t \equiv & 2s(u_{2} - 0.5) \\ z \equiv & s \frac{0.46341 t − 0.73369 t^{2} + 0.27385 t^{3}} {0.597999 − t + 0.309420 t^{2} + 0.093073 t^{3}} \\ y_{\tilde{m}} \equiv & z \sqrt{1 + x_{\tilde{m}}^{2}} \end{align} \]

ここで,\(u_{2}\)\([0, 1)\) の乱数である.

Rotate \(\omega_{m}\)

\[ \begin{align} \left( \begin{array}{c} x_{\tilde{m}} \\ y_{\tilde{m}} \end{array} \right) = \left( \begin{array}{c} cos \phi_{\omega_{i}'} & -sin \phi_{\omega_{i}'} \\ sin \phi_{\omega_{i}'} & cos \phi_{\omega_{i}'} \end{array} \right) \left( \begin{array}{c} x_{\tilde{m}} \\ y_{\tilde{m}} \end{array} \right) \end{align} \]

Unstretch \(\omega_{m}\)

\[ \begin{align} \left( \begin{array}{c} x_{\tilde{m}} \\ y_{\tilde{m}} \end{array} \right) = \left( \begin{array}{c} \alpha x_{\tilde{m}} \\ \alpha y_{\tilde{m}} \end{array} \right) \end{align} \]

Normalize \(\omega_{m}\)

\[ \begin{align} \omega_{m} \equiv \frac{(-x_{\tilde{m}}, -y_{\tilde{m}}, 1)} {\sqrt{x_{\tilde{m}}^{2} + y_{\tilde{m}}^{2} + 1}} \end{align} \]

まとめ

パラメータ

名称 記号
表面粗さ \(roughness\) \(0 < roughness \le 1\)

レンダリング

名称 記号
BRDF \(f_{r} \left( x, \omega_{i}, \omega_{o} \right)\) \(\frac{1}{4 \left( \omega_{i} \cdot \omega_{n} \right) \left( \omega_{o} \cdot \omega_{n} \right)} F \left( \omega_{i}, \omega_{m} \right) G_{2} \left( \omega_{i}, \omega_{o}, \omega_{m} \right) D \left( \omega_{m} \right)\)
pdf \(p \left( \omega_{o} \right)\) \(\frac{1}{4 \left( \omega_{i} \cdot \omega_{m} \right)} D_{\omega_{i}} \left( \omega_{m} \right)\)
反射レイのWeight \(\frac{f_{r} \left( x, \omega_{i}, \omega_{o} \right) \left( \omega_{o} \cdot \omega_{n} \right)} {p \left( \omega_{o} \right)}\) \(F \left( \omega_{i}, \omega_{m} \right) G_{1} \left( \omega_{o}, \omega_{m} \right)\)

以上,GGX BRDFについてまとめた.

参考文献

  1. B. Walter, S. R. Marschner, H. Li, K. E. Torrance: Microfacet models for refraction through rough surfaces (2007)
  2. E. Heitz, E. D'Eon: Importance Sampling Microfacet-Based BSDFs using the Distribution of Visible Normals (2014)
  3. B. Burley: Physically-Based Shading at Disney (2012)

2015年11月29日日曜日

Smooth Diffuse Surface (Lambert BRDF)

Smoth Diffuse Surface (Lambert BRDF)

本ページでは,モンテカルロレイトレーシング法において, 理想拡散反射表面をレンダリングする方法を説明する. 理想拡散反射表面として Lambert BRDF の実装に必要な パラメータや数式,重点的サンプリングの方法を説明する.

Lambert BRDF

Lambert BRDFは,全ての方向に放射輝度が等しく反射するBRDFで

\[ \begin{align} f_{r} \left( x, \omega_{i}, \omega_{o} \right) = \frac{\rho}{\pi} \tag{1} \end{align} \]

で定義される.ここで,\(\rho\) は表面の反射率を表す.

重点的サンプリング (Importance Sampling)

反射方向 \(\omega_{o}\) の重点的サンプリングを行う. Lambert BRDF では, \(\omega_{o}\) は反射エネルギー分布に従ってサンプリングする. すなわち,反射方向のpdf \(p \left( \omega_{o} \right)\)

\[ \begin{align} p \left( \omega_{o} \right) = \frac{cos \theta_{\omega_{o}}}{\pi} \tag{2} \end{align} \]

となる.ここで, \(cos \theta_{\omega_{o}} = \omega_{o} \cdot \omega_{n}\) である.

重点的サンプリングの際の反射レイのWeightは

\[ \begin{align} \frac{f_{r} \left( x, \omega_{i}, \omega_{o} \right) \left( \omega_{o} \cdot \omega_{n} \right)} {p \left( \omega_{o} \right)} = \rho \tag{3} \end{align} \]

となる.

\(\omega_{o}\) のサンプリング

反射方向 \(\omega_{o}\) のサンプリングには逆関数法を用いる.

注意として,このサンプリング方法は 法線が \(\omega_{n} = (0, 0, 1)\) の時の 反射方向をサンプリングしているため, サンプリング後は任意の法線に対応するように基底変換をする.

\[ \begin{align} P \left( \theta, \phi \right) = & \int_{0}^{\phi} \int_{0}^{\theta} p \left( \theta^{\prime}, \phi^{\prime} \right) \ sin \theta^{\prime} \ d \theta^{\prime} d \phi^{\prime} \\ = & \frac{\phi}{2 \pi} \cdot (1 - cos^{2} \theta) \\ P \left( \theta \right) = & 1 - cos^{2} \theta \\ P \left( \phi \right) = & \frac{\phi}{2 \pi} \end{align} \]

上記の \(P \left( \theta \right)\)\(P \left( \phi \right)\) から, 反射方向の \(\theta\)\(\phi\) は以下の式から求めることができる.

\[ \begin{align} \theta = & cos^{-1} (\sqrt{1 - u_{1}}) \tag{4} \\ \phi = & 2 \pi u_{2} \tag{5} \end{align} \]

ここで, \(u_{1}\)\(u_{2}\)\([0, 1)\) の乱数である.

まとめ

パラメータ

名称 記号
反射率 \(\rho\) \(0 \le \rho \le 1\)

レンダリング

名称 記号
BRDF \(f_{r} \left( x, \omega_{i}, \omega_{o} \right)\) \(\frac{\rho}{\pi}\)
pdf \(p \left( \omega_{o} \right)\) \(\frac{cos \theta_{\omega_{o}}}{\pi}\)
反射レイのWeight \(\frac{f_{r} \left( x, \omega_{i}, \omega_{o} \right) \left( \omega_{o} \cdot \omega_{n} \right)} {p \left( \omega_{o} \right)}\) \(\rho\)

以上,Lambert BRDFについてまとめた.

参考文献

  1. H. W. Jensen: Monte Carlo Ray Tracing (2003)
  2. H. Son: Realistic Image Synthesis with Light Transport (2015)

2015年11月23日月曜日

双方向反射分布関数 (Bidirectional Reflectance Distribution Function: BRDF)

双方向反射分布関数 (Bidirectional Reflectance Distribution Function: BRDF)

BRDFの定義についてメモする. BRDFの実装やそのテストを実装する際は, 実装内容がBRDFの定義を満たしているかを確認することで バグを発見しやすくなると思われる.

BRDFとは

  • 表面の見た目を決定づける役割を果たす
  • 入射光に対する反射光の方向の分布を表す

BRDFの定義

BRDFは数学的には放射照度に対する出射光の放射輝度の比を表し, 以下の関数 \(f_{r}\) で定義される.

\[ \begin{align} f_{r}(x, \boldsymbol{\omega_{o}} \leftarrow \boldsymbol{\omega_{i}}) \equiv & \frac{dL_{r}(x \rightarrow \boldsymbol{\omega_{o}})} {dE(x \leftarrow \boldsymbol{\omega_{i}})} \\ = & \frac{dL_{r}(x \rightarrow \boldsymbol{\omega_{o}})} {L_{i}(x \leftarrow \boldsymbol{\omega_{i}}) \ d \sigma_{p}(\boldsymbol{\omega_{i}})} \end{align} \]

ここで, \(\boldsymbol{\omega_{i}}\)\(\boldsymbol{\omega_{o}}\) は それぞれ光の入射方向・反射方向を表す単位ベクトルであり, 点 \(x\) での法線ベクトルを \(\boldsymbol{n_{x}}\) とするとき, \((\boldsymbol{n_{x}} \cdot \boldsymbol{\omega_{i}}) \geq 0\) を満たす方向ベクトル \(\boldsymbol{\omega_{i}} \in \Omega_{x}\) である( \(\boldsymbol{\omega_{o}}\) に関しても同様). \(L_{i}, L_{r}\) はそれぞれ点 \(x\) に入射, 点 \(x\) から反射する放射輝度を表す. \(\sigma_{p}\)\(\Omega_{x}\) の投影立体角測度であり, \(\Omega_{x}\) の立体角測度 \(\sigma\)

\[ \begin{align} d\sigma_{p}(\boldsymbol{u}) \equiv cos \theta_{\boldsymbol{u}} \ d\sigma(\boldsymbol{u}), \ \ \boldsymbol{u} \in \Omega_{s} \end{align} \]

を満たす関係である.ここで, \(cos \theta_{\boldsymbol{u}} = (\boldsymbol{n_{x}} \cdot \boldsymbol{u})\) である.

物理ベースなBRDF

物理ベースなBRDFは,反射する放射エネルギーが負になることはないので

\[ \begin{align} \forall \boldsymbol{\omega_{i}} \ \forall \boldsymbol{\omega_{o}}, \ f_{r}(x, \boldsymbol{\omega_{o}} \leftarrow \boldsymbol{\omega_{i}}) \ cos \theta_{\boldsymbol{\omega_{o}}} \ge 0 \end{align} \]

\(cos \theta_{\boldsymbol{\omega_{o}}} \ge 0\) であるため,

\[ \begin{align} \forall \boldsymbol{\omega_{i}} \ \forall \boldsymbol{\omega_{o}}, \ f_{r}(x, \boldsymbol{\omega_{o}} \leftarrow \boldsymbol{\omega_{i}}) \ge 0 \end{align} \]

である. また,物理ベースなBRDFは以下の2つの特性を持つ.

ヘルムホルツの相反性 (Helmholtz Reciprocity)

BRDFの値は入射・出射方向が入れ替わっても変化しない.

\[ \begin{align} f_{r}(x, \boldsymbol{\omega_{o}} \leftarrow \boldsymbol{\omega_{i}}) = f_{r}(x, \boldsymbol{\omega_{i}} \leftarrow \boldsymbol{\omega_{o}}) \end{align} \]

この特性があるため,Bidirectional Path TracingやPhoton Mappingのように 光源側とカメラ側の双方向から光の経路をトレースできる.

エネルギー保存則 (Energy Conservation)

ある点から反射した総出射光束は総入射光束以下となる.

\[ \begin{align} \rho(x) = \frac{\int_{\Omega} L(x \rightarrow \boldsymbol{\omega_{o}}) \ d \sigma_{p}(\boldsymbol{\omega_{o}})} {\int_{\Omega} L(x \leftarrow \boldsymbol{\omega_{i}}) \ d \sigma_{p}(\boldsymbol{\omega_{i}})} \le 1 \end{align} \]

この式から以下の式を導くことができる

\[ \begin{align} \forall \boldsymbol{\omega_{i}}, \ \int_{\Omega} f_{r}(x, \boldsymbol{\omega_{o}} \leftarrow \boldsymbol{\omega_{i}}) \ d \sigma_{p}(\boldsymbol{\omega_{o}}) \le 1 \end{align} \]

この式から,エネルギー保存則を満たすBRDFは 入射した放射エネルギー以上の放射エネルギーを反射しないことがわかる.

まとめ

BRDFの定義

\[ \begin{align} f_{r}(x, \boldsymbol{\omega_{o}} \leftarrow \boldsymbol{\omega_{i}}) \equiv \frac{dL_{r}(x \rightarrow \boldsymbol{\omega_{o}})} {L_{i}(x \leftarrow \boldsymbol{\omega_{i}}) \ d \sigma_{p}(\boldsymbol{\omega_{i}})} \end{align} \]

物理ベースなBRDF

名称
非負 \(\forall \boldsymbol{\omega_{i}} \ \forall \boldsymbol{\omega_{o}}, \ f_{r}(x, \boldsymbol{\omega_{o}} \leftarrow \boldsymbol{\omega_{i}}) \ge 0\)
ヘルムホルツの相反性 \(f_{r}(x, \boldsymbol{\omega_{o}} \leftarrow \boldsymbol{\omega_{i}}) = f_{r}(x, \boldsymbol{\omega_{i}} \leftarrow \boldsymbol{\omega_{o}})\)
エネルギー保存則 \(\forall \boldsymbol{\omega_{i}}, \ \int_{\Omega} f_{r}(x, \boldsymbol{\omega_{o}} \leftarrow \boldsymbol{\omega_{i}}) \ d \sigma_{p}(\boldsymbol{\omega_{o}}) \le 1\)

以上,BRDFの定義についてまとめた. BRDFの実装やテストの実装の際は, これらの定義を満たしているかを確認することでバグを発見しやすくなると思われる.

参考文献

  1. R. Montes, C. Ureña: An Overview of BRDF Models (2012)
  2. H. Son: Realistic Image Synthesis with Light Transport (2015)