Blenderから自作レンダラー向けに設定をエクスポート

こんにちは。 このページは レイトレ合宿6 の記事です．

今年もレイトレ合宿の日が近くなってきて どんなシーンをレンダリングしようかと考えている頃ではないでしょうか。 そこで、 ¹ Blenderで扱うシーンの設定を自作のレンダラー向けにエクスポートできるように してみます。 Blenderからエクスポートすることの利点は、

• Blenderという強力なGUIツールを用いて効率良くシーンを作成できる
• Cyclesレンダラーなどの他のエンジンと同じ設定でレンダリングすることで、レンダリング結果を比較できる

などが挙げられます。

デモファイルのダウンロード

設定をエクスポートするに当たって、 Blenderの公式サイトで配布されている ヘリコプターデモ “scene-Helicopter-27.Blend” を用います。

Blenderデモ

以下の画像の場所からダウンロードできます。感謝！。

Blenderの設定

設定をエクスポートする前に、Blenderの設定を確認します。

1. RenderEngineの設定

使用するレンダリングエンジンを選択します。 ここは Cycles Render を選択します (ヘリコプターデモの場合は既にCycles Renderが選択されています)。

2. Dimensionsの設定

出力する画像の解像度を設定します。 Resolution の項目を100%にし、XとYに画像の解像度を設定します。

3. Samplingの設定

Pattern で使用する乱数のサンプラーを選択し、 Settings の S に乱数のシード値を指定しします。 Samples の Render に1ピクセル当りのレイのサンプル数を指定します。

4. LightPathsの設定

様々な経路長のパスを考慮するため、Full Global Illumination を選択します。

5. PostProcessingの設定

画像処理前のレンダリング画像を取得するため、 Post Processingのチェックを全て外しておきます。

6. Surfaceの設定

Cyclesエンジンの初期設定で、うっすらとした環境光が設定されている場合があります。 意図しない光源が設定されていると比較の際に困るので、 Background の Strength が 0 になっているか確認しておきます。

7. Performanceの設定

パフォーマンスの比較を行う場合は、Threads で使用するスレッド数を指定します。

リファレンスの作成

以上の設定を確認してレンダリングを行います。

自作レンダラーとの比較のために、マテリアルの設定をもう少しシンプルにしておきます。

これで、設定をエクスポートする準備が整いました。

Blenderエクスポーターの作成

Blenderから設定をエクスポートするツールを作成します。

BlenderはPythonインタプリタを内蔵しており、 Blender内でPythonスクリプトを実行することができます。 そこで、Blenderから設定をエクスポートするPythonスクリプトを作成してみます。

Pythonスクリプトの実行

Pythonスクリプトは、Blender起動時に引数として渡して実行することができます。 試しにスクリプトを実行してみます。 以下は.Blendファイルのファイル名を表示するスクリプトです。

#
# exporter1.py
#

import bpy

file_name = bpy.path.display_name_from_filepath(bpy.data.filepath)
print("## FileName: {0}".format(file_name))

端末を開いて、Blenderに引数として渡して実行します。

% blender --background Helicopter.blend --python exporter1.py
Read blend: Helicopter.blend

## FileName: Helicopter

Blender quit

Blenderの実行ファイルの場所は、Windowsなら C:\Program Files\Blender Foundation\Blender\blender.exe、 macOSなら/Applications/Blender/blender.app/Contents/MacOS/blenderにあると思います。

--background オプションを指定すると、BlenderのGUIを起動せずに Pythonスクリプトを実行することができます。

実行する際のポイントは、Helicopter.blendよりも後に--python exporter1.py を指定することです。 Helicopter.blendの設定が読み込まれた後にスクリプトが実行されます。

スクリプト内でimport bpy を実行することで、BlenderのAPIを利用できるようになります。 利用できるAPIの詳細は Blender Documentation Python API で参照できます。

シーンデータの取得

.Blendファイルを読み込んだ際の各種設定は、 bpy.data に記録されます (bpy.dataの構造については BlendData を参照して下さい)。

各シーンの情報については bpy.data.scenes に記録されています。 (bpy.data.scenesの構造については Scene を参照して下さい)。

特定のシーンデータを取得する場合はシーン名を指定することで取得できます。 シーン名は Outliner で確認できます。

以下はSceneという名前のシーンから情報を取得するスクリプトの例です。

#
# exporter2.py
#

import bpy

# Get a scene
scene_name = "Scene"
scene = bpy.data.scenes[scene_name]

# Get a image resolution
x = scene.render.resolution_x
y = scene.render.resolution_y
print("Image resolution: {0} x {1}".format(x, y))

# Get a sampler type
sampler = scene.cycles.sampling_pattern
print("Sampler: {0}".format(sampler))

# Get a sampler seed
seed = scene.cycles.seed
print("Sampler seed: {0}".format(seed))

# Get the number of threads
threads = scene.render.threads
print("Threads: {0}".format(threads))

# Get a render samples
samples = scene.cycles.samples
print("Render samples: {0}".format(samples))

最初に、 scene = bpy.data.scenes[scene_name] で欲しいシーン情報を取得し、 シーンの各種情報を取得しています。 上のスクリプトでは、情報を標準出力に出力しているだけですが、 実際には自作レンダラーのフォーマットに加工してファイルに出力します。

実行結果は以下のようになります。

% blender --background Helicopter.blend --python exporter2.py
Read blend: Helicopter.blend

Image resolution: 1920 x 1080
Sampler: CORRELATED_MULTI_JITTER
Sampler seed: 123456789
Threads: 4
Render samples: 512

Blender quit

もし欲しい情報のAPIがわからない場合は、 BlenderのGUI上でその情報を編集するためのボタンやメニューの上に マウスカーソルを置いてみて下さい。 しばらくするとポップアップが現れて、その情報を利用するためのAPIを表示してくれます。

カメラ

シーンのアクティブカメラの情報は、scene.camera に入っています。

#
# exporter3.py
#

import bpy

# Get a scene
scene_name = "Scene"
scene = bpy.data.scenes[scene_name]

# Get an active camera
camera_obj = scene.camera

# Get rotation angles
for axis in camera_obj.rotation_euler.order:
    axis_index = ord(axis) - ord('X')
    angle = camera_obj.rotation_euler[axis_index]
    print("{0} axis rotation: {1} radian.".format(axis, angle))

# Get a camera position
location = camera_obj.location
print("position: ({0}, {1}, {2})".format(location[0], location[1], location[2]))

# Get a field of view
print("Horizontal field of view: {0}".format(camera_obj.data.angle_x))
print("Vertical field of view: {0}".format(camera_obj.data.angle_y))

scene.camera 自体は Object 構造になっており、 カメラの位置や回転などの情報が入っています。 scene.camera.data にはカメラの画角やレンズなどの情報が入っています (Camera を参照して下さい)。

実行結果は以下のようになります。

% blender --background Helicopter.blend --python exporter3.py
Read blend: Helicopter.blend

X axis rotation: 4.169113636016846 radian.
Y axis rotation: 3.1490345001220703 radian.
Z axis rotation: 3.9584245681762695 radian.
position: (1.483912467956543, -1.3075470924377441, 1.6609470844268799)
Horizontal field of view: 0.8575560450553894
Vertical field of view: 0.5033799409866333

Blender quit

オブジェクト

シーン内のオブジェクトのデータは scene.objects に入っています (各オブジェクトの構造は Object を参照して下さい)。 例えば、シーン内のオブジェクト名を全て表示するスクリプトは、

#
# exporter4.py
#

import bpy

# Get a scene
scene_name = "Scene"
scene = bpy.data.scenes[scene_name]

# Display object names
for obj in scene.objects:
  print(obj.name)

のようになります。

試しに上画像のOutlinerのシーンで実行してみると、

blender --background TestScene.blend --python exporter4.py
Read blend: TestScene.blend

ObjectB
ObjectA
ObjectA-2
ObjectA-1
Lamp
Camera

Blender quit

scene.objects にはオブジェクト同士の親子関係に関わらず全てのオブジェクトが 入っています。 親子関係を考慮する場合は、

#
# exporter5.py
#

import bpy

# Get a scene
scene_name = "Scene"
scene = bpy.data.scenes[scene_name]

def printChildObjects(objects):
  for obj in objects:
    print("  child node: " + obj.name)

# Display object names
for obj in scene.objects:
  if obj.parent == None:
    print(obj.name)
    printChildObjects(obj.children)

のように書くことができます。 obj.parent には親となるオブジェクトが入っており (トップレベルの場合は None)、 obj.children は子となるオブジェクトの配列(子がいない場合は要素を持たない配列) となっており、これらを用いて親子階層を意識した処理を書くことができます。

実行結果は

% blender --background TestScene.blend --python exporter5.py
Read blend: TestScene.blend
ObjectB
ObjectA
  child node: ObjectA-1
  child node: ObjectA-2
Lamp
Camera

Blender quit

となります。

次は、オブジェクト内の情報をエクスポートしてみます。

#
# exporter6.py
#

import bpy

# Get a scene
scene_name = "Scene"
scene = bpy.data.scenes[scene_name]

# Display object names
for obj in scene.objects:
  if not obj.type in {'MESH', 'CURVE', 'SURFACE'}:
    continue

  print(obj.name)

  # Get the active material of an object
  material = obj.active_material
  if material != None:
    print("  Material name: {0}".format(material.name))
    print("  Material index: {0}".format(bpy.data.materials.find(material.name)))

  # Export a object data to wavefromt format
  bpy.ops.object.select_all(action='DESELECT')
  obj.select = True
  bpy.ops.export_scene.obj( \
      filepath = obj.name + ".obj", \
      check_existing = True, \
      axis_forward = 'Y', \
      axis_up = 'Z', \
      use_selection = True, \
      use_animation = False, \
      use_mesh_modifiers = True, \
      use_edges = False, \
      use_smooth_groups = False, \
      use_smooth_groups_bitflags = False, \
      use_normals = True, \
      use_uvs = True, \
      use_materials = False, \
      use_triangles = True, \
      use_nurbs = False, \
      use_vertex_groups = False, \
      use_blen_objects = True, \
      group_by_object = False, \
      group_by_material = False, \
      keep_vertex_order = True, \
      global_scale = 1.0, \
      path_mode = 'AUTO')

obj.type を見ることで、どんな種類のオブジェクトなのかがわかります (例えば ‘MESH’、‘CURVE’、’CAMERA’など。Type of Object)。 また、この種類の違いによって、 obj.data の構造も異なります (’MESH’なら Mesh 、’CURVE’なら Curve など)。

bpy.ops.export_scene.obj はオブジェクトを .obj 形式で書き出すAPIです (Export Scene Operators を参照)。 オブジェクトの種類毎に処理を変えるのは大変なので、 ここでは bpy.ops.export_scene.obj を用いて ’MESH’や’CURVE’などのオブジェクトは全て .obj 形式に書き出すようにしました。

以下は実行結果の一部を切り出したものです。

% blender --background Helicopter.blend --python exporter6.py
Read blend: Helicopter.blend

...

Bolt+Nut.011
  Material name: Metal
  Material index: 3
    (  0.0409 sec |   0.0000 sec) OBJ Export path: 'Bolt+Nut.011.obj'
          (  0.2900 sec |   0.2481 sec) Finished writing geometry of 'Bolt+Nut.011'.
      (  0.2964 sec |   0.2550 sec) Finished exporting geometry, now exporting materials
      (  0.2964 sec |   0.2550 sec) OBJ Export Finished

Progress: 100.00%

...

Blender quit

オブジェクト“Bolt+Nut.011”が“Bolt+Nut.011.obj”として書き出されました。 自作レンダラーでは、このobjファイルを読み込むように設定をエクスポートしています。

マテリアル

マテリアルの情報は bpy.data.materials に入っています (各マテリアルの構造は Material を参照して下さい)。 Blenderのマテリアルデータは複雑で、また、自作レンダラーのマテリアルデータとは 取るパラメータも異なるため、 マテリアルデータに関してはエクスポートは行わず後から手動で設定するようにしました。

Cyclesレンダラーと自作レンダラーの比較

エクスポートした設定を自作レンダラーに読み込ませて実際にレンダリングしてみました。

以下はCyclesレンダラーとの比較になります。

上画像がCyclesレンダラーでレンダリングした画像で、 下画像が自作レンダラーでレンダリングした画像です。 カメラ設定やオブジェクトの場所は一致させられたと思います。

しかし、レンダラー間でのマテリアルが取るパラメータの違いや係数の違い、 色空間設定の違いから、 結果を完全に一致させるにはまだ工夫が必要だと感じました。 一応今のままでも簡単な比較には使えそうです。

私のレンダラー向けのエクスポーターは、 blender_scene_exporter.py にあります。 参考にどうぞ。

参考

1. Blender Documentation Python API: https://docs.blender.org/api/current/

---

1. この記事内で使用したBlenderのバージョンは 2.79b です。↩

レンダリング方程式の変数変換

このページは レイトレ Advent Calendar 2017 の記事である． 本ページでは，レンダリング方程式の変数変換について考える．

記号について

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

記号	説明
\(\omega_{n}\)	単位法線ベクトル
\(\omega_{i}\)	入射方向の単位ベクトル
\(\omega_{o}\)	反射方向の単位ベクトル
\(L_{o}\)	ある点からの反射光の放射輝度
\(L_{e}\)	ある点が自発光している場合の放射輝度
\(L_{i}\)	ある点への入射光の放射輝度
\(S^{2}\)	半球上の微小立体角の集合
\(M\)	空間中にある物体上の微小面積の集合
\(f_{r}\)	BRDF

微小立体角を積分する形のレンダリング方程式は式1で表される．

\[ \begin{align} L_{o} \left( x, \boldsymbol{\omega_{o}} \right) = L_{e} \left( x, \boldsymbol{\omega_{o}} \right) + \int_{S^{2}} f_{r} \left( x, \boldsymbol{\omega_{i}}, \boldsymbol{\omega_{o}} \right) L_{i} \left( x, \boldsymbol{\omega_{i}} \right) cos \theta_{\omega_{n}} d \sigma \left( \boldsymbol{\omega_{i}} \right) \tag{1} \end{align} \]

式1では， 方向 \(\omega_{i}\) を中心とした微小立体角を半球上で積分することによって 放射輝度の計算を行っている．

一方で，入射方向 \(\omega_{i}\) ，出射方向 \(\omega_{o}\) の代わりに， 下図のように空間中の3点を接続する形のレンダリング方程式を考えることもできる．

3点接続のレンダリング方程式では， 点 \(x^{''}\) を中心とした微小面積を 空間中の物体上で積分することで放射輝度の計算を行う．

微小面積の積分では，点 \(x\) から半球を通して見える範囲の微小面積のみを積分したい． そこで次の関数を導入する．

\[ \begin{align} V \left( x \leftrightarrow x^{''} \right) = \begin{cases} 1 \ \ \left( x と x^{''} の間に遮蔽物が無い場合 \right) \\ 0 \ \ \left( x と x^{''} の間に遮蔽物がある場合 \right) \end{cases} \tag{2} \end{align} \]

この関数によって，点 \(x\) から半球を通して見えない微小面積をカットする．

次に，微小立体角 \(d \sigma\) と微小面積 \(dA\) の関係は，下図から式3のように表すことができる．

\[ \begin{align} d \sigma = \frac{dA cos \theta_{\omega_{n^{''}}}}{\left| x - x^{''} \right|^{2}} \tag{3} \end{align} \]

式2，3から，式1を基にして微小面積を積分する形のレンダリング方程式は

\[ \begin{align} L_{o} \left( x \rightarrow x^{'} \right) = L_{e} \left( x \rightarrow x^{'} \right) + \int_{M} f_{r} \left( x^{'} \rightarrow x \rightarrow x^{''} \right) L \left( x \leftarrow x^{''} \right) V \left( x \leftrightarrow x^{''} \right) \frac{cos \theta_{\omega_{n}} cos \theta_{\omega_{n^{''}}}} {\left| x - x^{''} \right|^{2}} dA \left( x^{''} \right) \tag{4} \end{align} \]

となる．

まとめ

本ページではレンダリング方程式の変数変換について簡単にまとめた． 積分の変数変換ができれば，論文の式の検証ができたり自分で新しい式を考えることができて 便利です．

参考文献

1. James T. Kajiya: The rendering equation (1986)

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)

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の反射分布が正しいかというところまでは確認できていません． ここからテスト項目を増やしていって，確認できることを広げていこうと思います．

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)

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)

双方向反射分布関数 (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)