情報処理学会研究報告
IPSJ SIG Technical Report

メトロポリス光輸送法における双方向変異法のための
光路サンプリング戦略決定について
大津 久平1,2,a)

楽 詠コウ3侯 啓明4岩崎 慶5,2

土橋 宜典6,2

須田 礼二1

西田 友是2,7

概要：写実的な画像を生成するための大域照明計算はコンピュータグラフィックスにおいて重要な分野である. 近年, マルコフ連鎖モンテカルロ法を用いた大域照明計算法が活発に研究, 開発されている. 例えば, マルコフ連鎖モンテカルロ法 (MCMC) を用いるレンダリング手法の一つであるメトロポリス光輸送法(MLT) はある規則に従って光路を変異させながら，光路空間におけるエネルギーの寄与分布に従うようなマルコフ連鎖を生成する. MCMC を用いるレンダリング手法では, サンプルが次元の境界を超える際に次元間変異が用いられる. 例えば, MLT では次元間変異として双方向変異法を用い，光路の長さや構造を変異させる．このような次元間変異手法はしばしばサンプル過程のエルゴード性を保つために用いられ, い
くつかの手法にとっては次元間変異は必須である. しかしながら MCMC を用いたレンダリングにおいて,
次元間変異の効率を改善するための研究はほとんど試されてこなかった. そこで, 本稿では, 双方向変異法の効率の改良について考察する. 重点的サンプリングに類似した考え方に基づく双方向変異法の効率についての解析から, 光路サンプリング戦略の決定を行う分布の制御によって, 双方向変異法の採択確率を向上できることが分かる. この解析に基づき, 我々は焼きなまし法を用いて採択確率の期待値を最大化することで光路サンプリング戦略の決定分布を効率化する手法を提案する. この手法を用いることにより, 双方向変異法の特徴を残しながら高い遷移確率で光路をサンプリングすることができる.

1. はじめに間変異と呼ぶ. 例えば, MLT 法においては双方向変異法がそれにあたり, 光源と視点から双方向に光路を構築するこ

写実的な画像の生成はコンピュータグラフィックスにおとで反射回数を含む光路の構造を変化させる. このようないて重要な分野である. 写実性を実現するために, 相互反射

変異はサンプル過程のエルゴード性を保つために必須であを考慮する大域照明手法が特に重要である. 大域照明手法

るが, 次元を超えて変異するために光路空間において大域を用いることにより, レンダリング画像の質を向上させる的な変異を行う必要があり変異ステップにおいて採択確率ことができる. 近年, マルコフ連鎖モンテカルロ法 (MCMCが低く, しばしば低効率の要因となる. しかしながら, 現在法) を用いるレンダリング手法が盛んに研究, 開発されていまでに MCMC 法を用いるレンダリングにおいて次元間変る. 例えば, MCMC 法を用いるレンダリング手法のひとつ

異の効率を改善する手法は行われてこなかった. そこで, 本であるメトロポリス光輸送法 (MLT 法)[10] では, ある規則

稿では, 双方向変異法の効率の改良について考察する. 提案に従って光路を変異させることによってマルコフ連鎖を生法では重点的サンプリングに基づく考えによって光路サン成する. これによって光路のエネルギー分布に従うサンプ

プリングの分布を改良することで採択確率が改善できるこ

ル列を生成することができ, 独立なサンプル列を用いる場とを示し, 焼きなまし法を用いた最適化により採択確率の合よりもより効率的なレンダリングが可能となる. MCMC

推定値を最大化するような光路サンプリングの分布を求め法を用いるレンダリング手法において, 反射回数が変化す

る. 我々は提案手法を 2 次元の大域照明フレームワーク [4]るような変異を行うケースがある. このような変異を次元に実装し, いくつかの実験により提案法の効果を確認した.

1
2
3
4
5
6
7
a)

東京大学
UEI リサーチ
コロンビア大学
浙江大学
和歌山大学
北海道大学
広島修道大学h-otsu@is.s.u-tokyo.ac.jp

c 1959 Information Processing Society of Japan⃝2. 光路を用いた大域照明計算の定式化
本稿で用いる,Veach による大域照明計算の定式化 [8] について述べる. MLT 法はこの定式化によって説明される.
大域照明計算はそれぞれのピクセル j において観測されたエネルギー Ij を求めることによって定式化される.

1

情報処理学会研究報告
IPSJ SIG Technical Report

∫
Ij =

P

fj (x)dµ(x)(1)

次に MH 法について述べる. 現在のサンプルを Xi とするとき, まず試験的なサンプル Xi′ を遷移関数 T から生成

ここで P は光路空間であり, シーン表面 M ⊂ R3 上のする: Xi′ ∼ T (Xi → ·). T は遷移カーネル K とは異なる分点列によって定義される光路の集合として定義される. す布であることに注意する. このとき次のサンプル Xi+1 は

なわち,

P=

∞
∪

Pn

(2)

n=2

i

であり, ここで Pn = {x1 x2 . . . xn | x1 , . . . , xn ∈ M} はそれぞれの頂点数における光路の集合である. また µ は光路空間に対して定義される測度であり, ある D ⊂ P に対して,∫∑∞µ(D) = k=2 µk (D ∪ Pk ), µk (D) = D dA(x1 ) · · · dA(xk )と定義される. また fj は光路のエネルギーの寄与関数であり, 次のように定義される.

n−1∏with probability a(Xi → Xi′ )(6)otherwiseで定義される. ここで a(Xi → Xi′ ) を採択確率と呼び,

a(Xi → Xi′ ) = min(f (X ′ )T (X ′ → X)
1,f (X)T (X → X ′ ))(7)で定義される. MH 法の遷移カーネルは詳細つりあい条件
を満たし, f は定常分布となる.
収束性について, MH 法によって生成されるマルコフ連

fj (xn ) = Le (x1 → x2 )·Xi
Xi+1 =
X ′鎖が f が 0 ではない初期状態の選択に関わらず有限回の移G(xk−1 ↔ xk )fs (xk−1 → xk → xk+1 )動で状態空間 Ω の任意の場所へ到達できる (f -irreducible)

k=2· G(xn−1 ↔ xn )We(j) (xn−1 → xn )

ならば

(3)
1 ∑
M →∞
h(Xi ) −−−−→
a.s.
M i=1
M

ここで, fs は双方向反射率分布関数 (BSDF), Le は光源(j)の放射輝度関数で, We は j 番目のピクセルの応答関数∫
h(x)f (x)dµ(x)(8)Ω(センサーの強度) である. また G はジオメトリ項と呼ばが成り立ち, これは独立なケースで言うところの大数の強れ, 次のように定義される.

法則にあたる ([7] の Theorem 7.4 を参照).

G(x ↔ x′ ) = V (x ↔ x′ )|cos (θ0 ) cos (θ1 )|∥x − x′ ∥2(4)

ここで θ0 , θ1 はそれぞれ x, x′ における法線と x′ − x,x − x′ とのなす角である. V は可視関数であり, x, x′ がお互いに可視であれば V (x ↔ x′ ) = 1, 不可視であればV (x ↔ x′ ) = 0 で定義される.

MCMC 法をレンダリングに応用した初めての手法である.
MLT 法では MCMC 法として第 3 節で紹介した Metropolis空間は P であり, 光路の構造を直接変更することで変異を

本節では MLT 法の用いる, メトロポリス・ヘイスティン
グス法 (MH 法)[3], [6] について述べる. 導入として, まずは MCMC 法について説明を与える.

MCMC 法はある規則に従いサンプルを次々と生成していくことで, ある分布に従うサンプル列を生成する手法である. 状態空間を Ω とすると, MCMC 法は遷移カーネルと呼ばれる確率密度関数 K を用いて現在のサンプル Xi ∈ Ω から次のサンプル Xi+1 ∼ K(Xi → ·) を選択する. このように Xi+1 が Xi のみに依存する確率変数の列 X1 , X2 , · · · ∈ Ωをマルコフ連鎖と呼ぶ.
すべての i = 1, 2, . . . に対して Xi ∼ π のとき Xi+1 ∼ πならば, π は定常分布と呼ぶ. 遷移カーネル K に対して, す
べての x, y ∈ Ω に対して次の条件を満たすような π が存在するとき, K は詳細つりあい条件を満たすといい, このとき π は定常分布となる.

c 1959 Information Processing Society of Japan⃝

本節ではメトロポリス光輸送法 (MLT 法)[10] について述べる. MLT 法は Veach と Guibas によって開発された,

Hastings 法 (MH 法)[3], [6] を用いる. MH 法で用いる状態

3. メトロポリス・ヘイスティングス法

K(x → y)π(x) = K(y → x)π(y)4. メトロポリス光輸送法

行い, 寄与関数 f (式 3, j は省略) に従う光路からなるマル
コフ連鎖を生成する. MLT 法を Algorithm 1 に示す. ここで M は変異回数であり, f ilm はレンダリング結果を保存するヒストグラムである. なお f ilm は 3 次元のフレーム
ワークの場合 2 次元のヒストグラム (すなわち画像) であり, 2 次元のフレームワークの場合は 1 次元のヒストグラ
ムとなる.

4.1 光路の初期化
関数 InitialPath によって光路の初期化が行われる.

MCMC 法において, 初期状態の選択によって定常状態に至るまでにバイアスが存在することが知られており, これを start-up bias と呼ぶ. start-up bias を除去するために,
我々は Veach らの手法に倣って, 独立なサンプルを用いる
手法を用いて複数の光路を生成し, それらの寄与に応じた分布から初期サンプルを決定する. 今回独立な手法として(5)

双方向パストレーシング (BPT)[9] を用いた.

2

情報処理学会研究報告
IPSJ SIG Technical Report

図 1

MLT 法の変異手法. 光路の構造やシーン上の点の位置, 光線の角度を変更することに
よって変異を行う. 変異手法は反射回数が変化する双方向変異法と, 変化しない摂動法に大きく分けられる. 摂動法はさらに扱う光路の組み合わせによってより細かい手法に分かれる. ここで, L:光源, E:センサー, D:拡散反射面, S:鏡面である.

Algorithm 1 メトロポリス光輸送法
1: f ilm ← EmptyFilm()
2: x̄ ← InitialPath()
3: b ← EstimateNormalizingConstant()
4: for m = 0 to M − 1 do
5:
T ← SelectTransitionFunction()6:ȳ ← Sample( ȳ ∼ T (x̄ →)·)
7:(ȳ)T (ȳ→x̄)a ← min 1, ff (x̄)T(x̄→ȳ)8:
if Random() < a then9:x̄ ← ȳ10:
end if11:
AccumulateContribution(f ilm, x̄, M, b)
12:
m←m+1
13: end for
14: SaveFilm(f ilm)

図 2

独立な双方向変異法のステップ. 光路 x̄ の状態に関係なく光路 ȳ がサンプリングされる. 詳細は第 5 節を参照.

組み合わせを変えて変異するため光路空間において大域的な変異を担当する. このような変異はサンプル過程のエル

4.2 正規化定数の推定

ゴード性を保つために必須であるが, 変異ステップにおい寄与の計算に用いるための正規化定数 B :=目的の寄与関数 f は（一般に）正規化されていないため,∫
f (x̄)dµ(x̄)Ω

では双方向変異法の効率を改善すること着目する. 一方摂を推定する必要がある. ここでは正規化定数の推定値を動法 (perturbations) は反射の回数を変えない変異を担当て採択確率が低く, しばしば低効率の要因となる. 本研究

b := ⟨B⟩ とおく. なお, 採択確率 (式 7) の計算は f が正規

し, 双方向変異法と比較してよりローカルな変異を担当す化されているか否かは関係しないことに注意する. 我々はる. 光路の性質によって詳細な手法がいくつか提案されてBPT を用いて正規化定数の推定を行った.

いるが, 本稿では扱わない.

5. 独立な双方向変異法

4.3 寄与の算出
それぞれの変異のステップにおいて f ilm に足し合わせ

本研究では双方向変異法の大域的な変異により着目し,る寄与を計算する. x̄ ∈ P が定常状態にあると仮定するとかつ実験, 実装をシンプルにするためにオリジナルの双方

1 f (x̄)x̄ ∼ f /b であり, f ilm に与える寄与は M
f (x̄)/b

向変異法を簡略化した手法を導入し, これを用いて比較実

=

b
M

と

なる.

装, 実験を行った. 本節では新たに導入した双方向変異法である独立な双方向変異法について説明を与える. ここでの4.4 変異手法

独立とは遷移関数が現在の状態に依存しないことを意味す

MLT 法は様々な変異方法を組み合わせて実装される. そる. 独立な手法における遷移関数 Tbd に基づいて新たな光

れらは主に 2 種類に分類できる (図 1). ひとつは双方向変路 ȳ を生成する手順について説明する (図 2). まず, 現在の異法 (bidirectional mutation) であり, 光源と視点から双方

光路 x̄ の頂点をすべて消去する (図 2 (1)). 本来光路の頂向に光路を構築することで反射回数を含む光路の構造を変点の一部が選択され, 部分的に光路の頂点が消去されるが,化させることで変異を行う手法である. 反射回数や光路の独立な手法ではすべてを消去する. 次に, 新たに生成する

c 1959 Information Processing Society of Japan⃝3

情報処理学会研究報告
IPSJ SIG Technical Report

光路の頂点数 n′ を分布 pn := U{2, nmax } (範囲 [2, nmax ]の離散一様分布) からサンプリングする (図 2 (2)). ここで

nmax はユーザーが指定する光路の最大頂点数である. 本来の手法では pn として x̄ の頂点数に依存した分布を導入するが, 独立な手法では一様なものとした. 次に nmax に応じて定義される, 光源側, 視点側から生成する頂点の個数を振り分けるために, 光源側からの頂点数を選択する分布を psel (s; n′ ) := U{0, n′ } で定める. この分布は光路のサン

図 3

重点的サンプリングによる採択確率の向上. 左図:分布がよりプリング戦略を決定する分布となる. これにより光源側か平らになることにより, 採択確率を向上できる. 図中では青線ら生成する頂点数 s ∼ psel (s; n′ ), 視点側から生成する頂点の R よりも赤線の R のほうが望ましい. 右図: f が 0 の領域′

数 t := n − s が定まる (図 2 (3)). 最後に, 光源側, 視点側から延びる部分光路をそれぞれ頂点数 s, t だけサンプリン
グして, 端点を繋ぐことによって新たな光路 ȳ を生成する(図 2 (4)). このとき, 遷移関数 Tbd は次のようになる.
∑
Tbd (x̄) = pn (nx̄ )psel (s; nx̄ )ps,t (x̄) (9)1=nmax nx̄ + 1の関係については [1], p.8 等を参照). しかしながら, 提案法

ps,t (x̄)

法ではなく, 大域的な変異を許容しながら採択確率の向上(10)を行う手法であり採択確率の向上による効果が期待できる.

s,t∈Z+ ,s+t=nx̄

= (x1 , . . . , xs ), x̄E=E(xn , . . . , xn−t ) とすると ps,t (x̄) := pL
s (x̄L )pt (x̄E ). ただE
し pL
s , pt はそれぞれ光源側, 視点側の部分光路をサンプリ
L
ングする確率であり pL
0 (x̄L ) = 1, p1 (x̄L ) = pA (x1 ). またs ≥ 2 のとき,·pσ⊥ (xi → xi+1 )G(xi ↔ xi+1 )

(11)
また

pσ⊥ (xn−i+1 → xn−i )G(xn−i+1 ↔ xn−i )

(
)
R(ȳ)
a(x̄ → ȳ) = min 1,
R(x̄)

(13)

を示す. 一般に重点的サンプリングによって R はより平た
布を取ることができる場合, R は定数関数となり, 採択確率は常に 1 となる. 例えば, 重点的サンプリングによって得られる R を R′ として, このとき光路空間上の 2 点 x̄, ȳ ∈ P

pEs (x̄E ) = pA (xn )
t−1
∏

に書ける.

いな分布となる (図 3 左). 例えば T として f に比例する分

t ≥ 2 のとき,·

f (x̄)

関数 R(x̄) := T (x̄) を用いて採択確率 (式 7) は次のよう

f に似た形とすることが採択確率の向上に有効であること

i=1
E
同様に視点側も, pE
0 (x̄E ) = 1, p1 (x̄E ) = pA (xn ).

6.1 重点的サンプリングによる採択確率の向上

R に対する重点的サンプリング, すなわち T をできるだけpLs (x̄L ) = pA (x1 )
s−1
∏

しい.

はカーネルサイズを変更することで採択確率を向上する方∑

こ こ で x̄ = (x1 , . . . , xn ), x̄Lできる. 図中では赤線の採択確率よりも青線の採択確率が望まくなり, 分散が増大する. (MCMC 法における相関と分散s,t∈Z+ ,s+t=nx̄

1のサンプリングを減らすことにより, 平均的に採択確率を向上を考え, R(x̄) > R(ȳ) > 0, R′ (x̄) > R′ (ȳ) > 0 であると

(12)

i=1

6. 双方向変異法の効率改善
本節では第 5 節で導入した, 独立な双方向変異法の効率を改善するための新たな手法を提案する. 双方向変異法の特徴として, 大域的な変異を担当するため変異ステップに仮定する. このとき, 重点的サンプリングによって R の値の相対差分 ∆R :=R(x̄)−R(ȳ)
がより小さくなる, すなわちR(x̄)

∆R > ∆R′ と仮定すると,
R′ (ȳ)R(ȳ)− ′= (1 − ∆R ) − (1 − ∆R′ )R(x̄) R (x̄)= ∆ R ′ − ∆R < 0

(14)
(15)

おいて採択確率が低く, しばしば低効率の要因となる. こであるから, 採択確率の定義における R の比率 R(ȳ)/R(x̄)の問題を緩和するため, 提案法では重点的サンプリングにはより大きくなる. また重点的サンプリングによって f の基づく考えによって光路サンプリングの分布を改良するこ

値が 0 の領域 (P \ supp (f )) をサンプリングする確率を減とで採択確率が改善できることを示し, 焼きなまし法を用らすことができれば, 平均的な採択確率が増加する (図 3

いた最適化により採択確率の推定値を最大化するような光

左). f の値が 0 の領域は採択確率が常に 0 であり, その領路サンプリングの分布を求める. 採択確率が大きいことは域に変異した光路は必ず棄却される.

必ずしもサンプリング効率の向上につながるわけではない.

さて, 独立な双方向変異法の Rbd は次のようになる. な

例えば, 摂動法においてはカーネルサイズを小さくするこお仮定より, pn は採択確率の計算の際にオフセットされるとで採択確率を向上できるが, サンプル同士の相関が大きので含めていないことに注意する.

c 1959 Information Processing Society of Japan⃝4

情報処理学会研究報告
IPSJ SIG Technical Reportf (x̄)
p(s; nx̄ )ps,t (x̄)
sels+t=nx̄Rbd (x̄) = ∑(16)

決定の最適化

双方向変異法の遷移確率は, 様々なサンプル戦略の分布を混合したものとみなすことができる. 本来の双方向変異法
において psel は一様分布であるが, 本手法では上での考察を踏まえて, 採択確率を分布 psel の選択を変更することで遷移確率 Tbd を f によりフィットしたものとして選択する
ことを試みる.

6.2 光路サンプリング戦略決定の最適化
本節では psel を決定するための最適化手法を提案する.
本稿では焼きなまし法 [2], [5] に基づいた最適化手法を用いる. 目的関数として以下で示される採択確率の期待値

a = E[a(x̄ → ȳ)] の推定値 ⟨a⟩ を用いた.
b
1 ∑a(x̄i → x̄i+1 )Mb i=1M⟨a⟩ =Algorithm 2 焼きなまし法による光路サンプリング戦略(17)

1: function OptimizeSelectionSA
2:
T ← Tinit
3:w←1
4:⟨a⟩ ← EstimateExpectedA(w)5:for i = 0 to Mopt − 1 do
6:
w′ ← Perturb w
7:⟨a′ ⟩ ← EstimateExpectedA(w′ )8:∆ ⟨a⟩ ← ⟨a′ ⟩ − ⟨a⟩9:if ∆ ⟨a⟩ ≥ 0 then10:w ← w′
11:⟨a⟩ ← ⟨a′ ⟩
12:
else
13:
if Random() < e−∆⟨a⟩/T then14:w ← w′
15:⟨a⟩ ← ⟨a′ ⟩16:
end if17:
end if18:
T ← T · (1 − rcool )19:
end for
20: end function

ここで x̄1 , x̄2 , . . . , x̄Mb は f に従う定常状態のマルコフ連鎖であり, Mb は最適化 1 ステップあたりの変異回数

Algorithm 3 採択確率の期待値推定である. 最適化ステップの疑似コードを Algorithm 2, 3

1: function EstimateExpectedA(w)2:x̄ ← InitialPath()3:
Aexp ← 0
4:
for i = 0 to Mb − 1 do
5:ȳ ← Sample( ȳ ∼ Tbd (x̄) → ·)に示す. ここで, Tinit は初期温度, rcool は冷却率, Mopt
は最適化のイテレーション回数である. また, wn,i :=psel (i; n) (i = 1, . . . , n), wn := (wn,1 , wn,2 , . . . , wn,n ),w := (w2 , w3 , . . . , wnmax ) として, 頂点数ごとに定義される psel をパラメータ化する. 我々はこの最適化ステップ
をレンダリング前に一度実行し, パラメータ w を決定してからレンダリングに移行する.

Rbd (ȳ;w)
6:a ← min 1, R
bd (x̄;w)
7:
Aexp ← Aexp + Mab
8:
end for9:
return Aexp
10: end function

6.2.1 採択確率の期待値推定
採択確率の期待値推定のステップ (Algorithm 3) について, パラメータ w を変更した場合, マルコフ連鎖をそのまま継続すると詳細つりあい条件を満たさず, 定常分布に従わなくなるため毎度呼ばれる毎に初期光路を再生成する必要がある. また今回目的関数が推定値のため, 最適化の結果得られた結果もまた推定値となる. 推定値を用いた最適化は実験的には正しく動作している, すなわち, 採択確率の期待値の推定値を最大化するようなパラメータが得られているが, 詳しい解析については今回の研究では行っておらず, 効率性や正しさについての議論は今後の課題とする.

6.2.2 パラメータ w の初期化と摂動
w の初期値はすべてのパラメータが 1 であるとした (Algorithm 2, 2 行目). パラメータがすべて 1 のとき, psel は
既存の双方向変異法と同等となる. また, 焼きなまし法においてパラメータを摂動させる必要がある (Algorithm 2,6 行目). パラメータ w のそれぞれの要素 w を次のように更新する. ただし, U は [0, 1] に分布する一様乱数である.

w ← max (0, min (1, w + (2U − 1) · T ))

(18)

7. 結果
7.1 実験環境
我々は提案法を 2 次元の大域照明レンダラに実装した.
これは Jarosz らの 2 次元における大域照明計算の定式化 [4]に基づいた実装である. 2 次元のフレームワークではシー
ンは 2 次元空間上に定義され, レンダリング結果は 1 次元のヒストグラムとして得られる. リファレンスとして双方
向パストレーシング (BPT)[9] を, また比較には独立な双方向変異法のみを用いる MLT 法を実装した. これらの実装
はすべて C++言語で行い, OpenMP による並列化を行った. またすべての実験は Intel i7-4930K CPU 3.4GHz のマ
シンで 12 スレッドを用いて行った.
実験に用いたシーンは簡単のため, 光源とカメラはエリ
アライト, エリアセンサーのみとしている. カメラはピン
ホールカメラではなく, 長さを持ったセンサーを直接シーン中に置いていることに注意する. これは 3 次元のフレームワークにおいてフィルムを直接シーン中に置くのと同等である. 図中の L は光源を表し, E はセンサーを表す. また, D は拡散反射面を表す. 図 4 (Scene1) はもっとも単c 1959 Information Processing Society of Japan⃝5

情報処理学会研究報告
IPSJ SIG Technical Report
0.95
Reference
MLT (uniform)
MLT (optimized)
0.9

0.85

0.8

図 4 Scene1: 2 本の線分から成るシーン. 2 本の線分は同一長かつ0.75

平行であり, 線分間の距離は線分の長さと等しい.
0.7

0.65

図 7

Scene1 を用いた推定された分布の比較. 横軸はセンサーの位
置, 縦軸は推定されたエネルギーである.
1.3

図 5

Scene2: 3 本の線分から成るシーン. 構成は Scene1 とほぼ

1.25

同様だが, 拡散反射面がひとつ追加されている.

Reference
MLT (uniform)
MLT (optimized)

1.2

1.15

1.1

1.05

1

0.95

0.9

図6

Scene3: 閉じた環境を持つシーン. シーン中心部にあるブロッカーによって, 光路は右回りか左回りかの経路で相互反射を行

0.85

図 8

う必要があり, 光路空間に複雑さを与えている.

Scene2 を用いた推定された分布の比較. 横軸, 縦軸は図 7 と同様.

純なシーンであり, 上部の光源 (Le = 1) と下部のセンサー

0.3
Reference
MLT (uniform)
MLT (optimized)(We = 1) だけで構成される. このシーンでは相互反射は発生せず, 光路空間の頂点数が 2 のものだけから成る. 図 5

0.28(Scene2) は Scene1 に拡散反射面 (拡散反射率 R = 1) となる線分を加えたもので, 加えた面を介して相互反射が発生するようなシーンである. 図 6 (Scene3) は閉じた環境の

0.26

0.24

シーンである. 中心に障害物を置き, L から E に到達する
ためには相互反射を介する必要がある. また光路が存在で0.22きる空間が左右に分かれているのも特徴である. このシー
0.2

ン中の表面は光源やセンサーを含め拡散反射率 R = 0.5 の拡散反射面である. すなわち, すべての面で拡散反射が起
こり得る.

0.18

図 9

Scene3 を用いた推定された分布の比較. 横軸, 縦軸は図 7 と同様.

7.2 実験
実験では提案法 (最適化したパラメータを用いた独立なストグラムに分割してエネルギー分布の推定を行った. 3

双方向変異法), 従来法 (独立な双方向変異法) に対して採次元のフレームワークではひとつのビンが 1 ピクセルに相択確率の期待値の推定と推定された分布（レンダリング結当する. 表 1 がそれぞれのシーンに対する推定された採択

果）の可視化を行い, 比較検討を行った. シーンごとに使確率の期待値である. この表からすべてのシーンにおいて用したパラメータを表 1 に示す. 提案法と従来法はレンダ

採択確率の期待値が上昇していることが確認できる. しか

リングに用いるサンプル数, シードサンプル数, 最大パス長

しかがら, シーンの複雑性が上がるにつれ上昇度が下がっの条件を揃えて実験を行った. センサーはビン数 100 のヒていることも確認できる. 図 7, 8, 9 がそれぞれのシーン

c 1959 Information Processing Society of Japan⃝6

情報処理学会研究報告
IPSJ SIG Technical Report
表 1

それぞれのシーンに対するパラメータ設定. Render はレンダリングのためのパラメー
ターで, Opt. は最適化のためのパラメーターである.

Scene

Render #mut

Render #seed sample

nmax

Opt. #iter (Mopt )Opt. #mut (Mb )Scene1

106

105

2

1000

50000

Scene2

106

105

3

5000

10000

Scene3

6

5

4

5000

10000

10

10
表 2

それぞれのシーンに対する推定された遷移確率の期待値.

Scene

Orignal E[a] (既存法)

Optimized E[a] (提案法)

Diﬀ.

Scene10.573250

0.889917

0.316667

Scene20.353590

0.539980

0.186390

Scene30.097555

0.129551

0.031996に対して推定された分布である. それぞれの図の横軸がセンサーの位置を表し, 縦軸が推定されたエネルギーを表す.
またリファレンス (緑線) は BPT を用いて 108 サンプルを用いてレンダリングした結果である.

[7]
[8]

8. まとめと今後の課題[9]

本稿では簡略化した場合の MLT 法の双方向変異法における効率を改善するための考察を行った. 提案法では重点的サンプリングに基づいた双方向変異法のにおける採択確[10]

lations by fast computing machines, J. Chem. Phys.,Vol. 21, No. 6 (1953).
Robert, C. P. and Casella, G.: Monte Carlo statisticalmethods, Springer-Verlag (2004).
Veach, E.: Robust Monte Carlo Methods for Light
Transport Simulation, Ph.D. thesis, Stanford university, Department of Computer Science (1997).
Veach, E. and Guibas, L.: Bidirectional estimator for
light transport, Computer Graphics (Proc. of SIGGRAPH 95), pp. 419–428 (1995).
Veach, E. and Guibas, L. J.: Metropolis light transport,
Proc. of SIGGRAPH 97, pp. 65–76 (1997).

率を向上させるための考察を行い, 焼きなまし法を用いた最適化によって光路サンプリング戦略の決定分布を選択することで, 双方向変異法の大域的な変異を保ちながら, 採択
確率を向上させる手法を提案した. 本稿では 2 次元の大域
照明フレームワークを利用した提案法の検証を行い, いく
つかのシーンに対して実際に採択確率が向上するような光路サンプリング戦略の決定分布を見つけ出した.
本稿では, 評価を簡単にするために, MLT の双方向変異法を簡略化したものを考えたが, より精密な議論のため, もともとの双方向変異法を用いた場合の評価もしていきたい.
参考文献
[1]

[2]

[3]

[4]

[5]

[6]

Brooks, S., Gelman, A., Jones, G. L. and Meng, X.-L.:
Handbook of Markov chain Monte Carlo, CRC Press(2011).
Cerny, V.: Thermodynamical approach to the traveling salesman problem: An eﬃcient simulation algorithm, Journal of Optimization Theory and Applications, Vol. 45, pp. 41–51 (1985).
Hastings, W. K.: Monte carlo sampling methods using markov chains and their applications, Biometrika,
Vol. 57, No. 1, pp. 97–109 (1970).
Jarosz, W., Schonefeld, V. and Kobbelt, L.: Theory,
analysis and application of 2D global illumination, ACM
Transactions on Graphics, Vol. 31, No. 5, p. Article 125(2012).
Kirkpatrick, S., Gelatt Jr., C. D. and Vecchi, M. P.: Optimization by Simulated Annealing, Science, Vol. 220,
No. 4598, pp. 671–680 (1983).
Metropolis, N., Rosenbluth, A. W., Rosenbluth, M. N.,
Teller, A. H. and Teller, E.: Equation of state calcu-

c 1959 Information Processing Society of Japan⃝

7