CG
光学
BRDF

はじめに

この記事はレンダリング方程式の BRDF について,特に鏡面BRDFについて説明します.最初はレンダリング方程式の定義から見直します.そこから Torrance-Sparrow モデル,Cook-Torrance モデル,Walter モデルと見ていきます.次に,Phong モデル,Blinn-Phong モデル,そして,異方性反射モデルとして,Ward モデル,Ashikhmin モデル,Schlick モデル,Disney モデルを取り上げます.最後に,微小面の法線分布とマスキング・シャドウイングについての補足をします.

BRDF

光は物体表面に当たると,散乱か吸収,もしくは両方の現象が発生します.散乱した光はまた別の物体表面に当たって...というのを繰り返して私たちの眼に入射してきます.光を光線と考え,ある物体表面上の点 $x'$ に光線が当たり,散乱した光線がまた別の物体表面上の点 $x$ に当たる,つまり,ある点 $x'$ から別の点 $x$ に光が伝達されていると考えることができます.また,ある点 $x'$ はさらに他の点 $x''$ から光が伝達されてきています.光が物体表面上で当たった時に,反射される放射エネルギーと入射したときの放射エネルギーとの比を反射率とし,これを計算することで光伝達のシミュレーションを行うことができます.レンダリング方程式ではこの反射率を表す関数として BRDF があります.BRDFは双方向反射率関数のことで,物体表面上に入射してきた放射エネルギーが,ある方向に出射する放射エネルギーの割合を表す関数です.BRDF には重要な2つの特性があります.1つ目がヘルムホルツの法則で,入射方向と出射方向を逆にしても,反射率が変わらないことです.これが双方向関数であることの理由です.2つ目がエネルギー保存則です.反射する放射エネルギー量は入射した放射エネルギー量を超えません.この2つの性質により,物理的な法則に基づいて計算することができます.BRDF は一般的に拡散BRDFと鏡面BRDFがありますが,今回は鏡面BRDFについて解説します.

レンダリング方程式

Kajiya[1986]によって,熱伝導の考え方をコンピュータグラフィックスに導入し,レンダリング方程式に統一しました.

I(x,x') = g(x,x')\left[
  e(x,x')+
  \int_{S}\rho(x,x',x'') I(x',x'') dx''
  \right]. \tag{1.1}

$S$ はシーン内のすべての表面上の点です.

$I(x,x')$ は表面の点 $x'$ から 点 $x$ に入る光の強さで,これを2点伝達強度といい,以降は伝達強度とします.伝達強度 $I$ は $x$ の単位面積,$x'$ の単位面積,単位時間時間当たりの放射エネルギーです.

dE = I(x,x') dtdxdx'. \tag{1.2}

$g(x,x')$ は幾何項で,表面の点が別の表面によって隠蔽されているかを表しています.$x'$ と $x$ が互いに見えない場合,幾何項はゼロになります.また,互いに見えるときの幾何項は $\frac{1}{r^2}$ です.$r$ は $x'$ から $x$ までの距離で $r = |x'-x|$ です.

$e(x,x')$ は放射光で表面の点 $x'$ から $x$ に出射される放射エネルギーです.これも $x$ の単位面積,$x'$ の単位面積,単位時間時間当たりの放射エネルギーなので

dE = \frac{1}{r^2}e(x,x') dtdxdx'. \tag{1.3}

$\rho(x,x',x'')$ は散乱項で,表面の点 $x''$ から $x'$ で散乱し,$x$ に入射する伝達強度の割合です.これを3点伝達反射率といい,以降は伝達反射率とします.この伝達反射率に点 $x''$ から $x'$ の伝達強度 $I(x',x'')$ を乗算します.

dE = \frac{1}{r^2}\rho(x,x',x'') I(x',x'') dtdxdx'dx''. \tag{1.4}

次にこれらの値と放射輝度やBRDFに結びつけます.
放射輝度は,単位時間・単位立体角・単位投影面積当たりのエネルギーです.点 $x$ に $(\theta',\phi')$ 方向から入射する放射輝度は

dE = i(\theta',\phi') d\omega dx'_{p} dt. \tag{1.5}

tuto-specular-brdf-two-point-transport.png

単位立体角 $d\omega$ は,微小面積 $dx$ を半径 $r$ の球上に投影した面積 $dx\cdot \cos\theta$ から

d\omega = \frac{1}{r^2}\cos\theta dx. \tag{1.6}

単位投影面積 $dx'_{p}$ は $dx'\cdot\cos\theta'$ だから,これらを式(1.5)に代入すると

dE = i(\theta',\phi') \frac{1}{r^2}\cos\theta\cos\theta' dtdxdx'. \tag{1.7}

式(1.2)と(1.5)が等しいとすると,伝達強度と放射輝度の関係は次のようになります.

I(x,x') = i(\theta',\phi')\frac{1}{r^2}\cos\theta\cos\theta'. \tag{1.8}

放射項も同じ関係で,

e(x,x') = e(\theta',\phi')\cos\theta\cos\theta'. \tag{1.9}

そして,伝達反射率を双方向反射関数(BRDF)に対応させます.BRDF $\rho(\theta',\phi',\psi',\sigma')$ は,ある点 $x'$ に入射してきた放射輝度を微小面積 $dx'$ に投影した放射エネルギー(これを放射照度ともいいます)と,そこから点 $x$ 方向に出射する放射輝度との比です.

i(\theta',\phi') = \rho(\theta',\phi',\psi',\sigma') i(\psi',\sigma') d\omega'' \cos\psi'. \tag{1.10}

tuto-specular-brdf-three-point-transport.png

立体角 $d\omega''$ は

d\omega'' = \frac{dx_{p}''}{r''^2} = \frac{1}{r''^2} \cos\psi''dx''. \tag{1.11}

式(1.11)の立体角 $d\omega''$ を式(1.10)に代入すると

i(\theta',\phi') = \rho(\theta',\phi',\psi',\sigma') i(\psi',\sigma') \frac{1}{r''^2}\cos\psi'\cos\psi''. \tag{1.12}

ここで

I(\psi',\sigma') = i(\psi',\sigma') \frac{1}{r''^2}\cos\psi'\cos\psi'', \tag{1.13}

とすると

i(\theta',\phi') = \rho(\theta',\phi',\psi',\sigma') I(\psi',\sigma'). \tag{1.14}

移項して BRDF を求める形にします.

\rho(\theta',\phi',\psi',\sigma') = \frac{i(\theta',\phi')}{ I(\psi',\sigma')}. \tag{1.15}

式(1.4)から $\rho(x,x',x'')$ を求めると

\rho(x,x',x'') = \frac{dE\cdot r^2}{I(x',x'')}. \tag{1.16}

今,

I(x',x'') = I(\psi',\sigma'), \tag{1.17}

なので,

\rho(x,x',x'') = \frac{dE\cdot r^2}{I(\psi',\sigma')}. \tag{1.18}

式(1.8)を代入すると

\begin{align}
\rho(x,x',x'') &= \frac{i(\theta',\phi') \frac{1}{r^2}\cos\theta\cos\theta' \cdot r^2}{I(\psi',\sigma')}\\[2ex]
&= \frac{i(\theta',\phi')}{I(\psi',\sigma')} \cos\theta\cos\theta'
\end{align} \tag{1.19}

式(1.15) から

\rho(x,x',x'') = \rho(\theta',\phi',\psi',\sigma') \cos\theta\cos\theta'
\tag{1.20}

点 $x$ に向かって反射される放射輝度(点 $x$ に入射していない状態)を考えれば $\cos\theta\cos\theta'$ は無視できます.また,点と点との間に遮蔽物がないとすれば,幾何項 $g$ もいらなくなります.点 $x'$ から 点 $x$ に反射する放射輝度は,点 $x'$ に半球上のあらゆる方向から入射してきた放射照度と BRDF を乗算したものを積分したものと,点 $x'$ から自己放射される放射輝度の和となります.これを表した式がレンダリング方程式です.

L_{o}(x,\vec{\omega}) = L_{e}(x,\vec{\omega}) + \int_{\Omega}
f_{r}(x,\vec{\omega}',\vec{\omega})
L_{i}(x,\vec{\omega}')
(\vec{\omega}'\cdot \vec{n}) d\vec{\omega}' \tag{1.21}

微小面反射モデル

Torrance-Sparrow モデル

物体表面に凹凸がある鏡面反射強度のピークは,正反射方向からグレージング角に向かって少しずれています.この性質をオフスペキュラーといいます.Torrance-Sparrow モデルは,物体表面が微小面の集合で構成されていると仮定し,微小面による隠蔽やフレネル反射を解析することでオフスペキュラーを表現できる物理ベースモデルです.

表面は等方性と仮定します.入射方向は表面の法線からの天頂角 $\psi$,反射方向は表面法線からの天頂角 $\theta$ と方位角 $\phi$ で表します.反射放射束のパラメータは入射方向と反射方向の3つの角度で表すことができます.

tuto-specular-brdf-torrance-sparrow-radiation.png

$dA$ に入射してくる放射輝度 $L_{i}$ は

L_{i} = \frac{d\Phi_{i}(\phi)}{\cos\psi dA d\omega_{i}}.

放射輝度 $L_{i}$ を使って,$dA$ に入射してくる放射束 $d\phi_{i}$ を表すと

d\Phi_{i}(\phi) = L_i(\psi)\cos\psi dA d\omega_{i}. \tag{2.1a}

$(\theta,\phi)$ 方向からの放射輝度は $dL_{r}(\psi;\theta,\phi)$.

dL_{r}(\psi;\theta,\phi) = \frac{
  d\phi_{r}(\psi;\theta,\phi)
}{
  \cos\theta dA d\omega_{r}
}. \tag{2.1b}

双方向反射率 $\rho$ は $(\theta,\phi)$ 方向からの反射放射輝度 ($dL_{r}$) から単位面積当たりの入射放射束($d\Phi_i/dA$) で割ったものです.

\rho(\psi;\theta,\phi) = \frac{dL_{r}(\psi;\theta,\phi)}{d\phi_{i}(\psi)/dA} = \frac{dL_{r}(\psi;\theta,\phi)}{L_{i}(\psi)\cos\psi d\omega_{i}}. \tag{2.2}

Torrance-Sparrow モデルは,表面の反射を拡散反射と鏡面反射の2つの成分で構成されています.

dL_{r}(\psi,\theta,\phi) = dL_{r,s}(\psi;\theta,\phi) + dL_{r,d}(\psi). \tag{2.3}

放射輝度 $N_{i}$ の光源が入射角 $\psi$ に向いているとき,拡散反射する光は表面の単位面積上に角度 $\psi$ で入射するので,$\cos(\psi)$ で変化します.なので

dL_{r,d}(\psi) = aL_{i}\cos\psi. \tag{2.4}

$a$ は任意の定数です.

下図に示す鏡のような単一な微小面からの放射束を考えます.放射束は微小面の法線 $n'$ に対して角度 $\psi'$ で表面に当たります.鏡面反射方向は方向 $\theta' = \psi'$ になります.通常は微小面上での方向 $\psi', \theta'$ は,表面上での方向 $\psi, \theta$ と一致しません.後者は単位面積 $dA$ の法線 $n$ が基準になっています.これらの角度を関連づけます.

tuto-specular-brdf-torrance-sparrow-microfacet.png

最初の図に示されているように,点 $O$ に $d\omega_{i}$ 方向から入射してくる放射輝度 $L_{i}$ の微小面積 $dA$ を考えます.$O$ における微小面の法線は入射ベクトル $d\omega_{i}$ と反射ベクトル $d\omega_{r}$ の中間ベクトル $d\omega'$ で,角度 $VOZ = \alpha$ だけ傾いています.角度 $QOV$ は微小面に対する入射角 $\psi'$ で,角度 $VOT$ はその微小面の反射角 $\theta'$ です.

$Pd\omega'$ は法線が $d\omega'$ に含まれる単位面積当たりの微小面の数とします.等方性表面の場合、表面法線 $OZ$ についての回転対称である $P$ のガウス分布と仮定します.

P = P(\alpha) = b \cdot e^{-c^2\alpha^2}, \tag{2.5}

ここで,$b$ と $c$ は定数です.立体角 $d\omega'$ の内側を向いている $dA$ 内の微小面の数は、

P(\alpha) d\omega' dA.

各微小面の面積を $f$ であるとすると、反射する微小面の全面積は、

fP(\alpha) d\omega' dA.

この領域の入射方向への投影面積は、

f \cos\psi' P(\alpha) d\omega' dA.

入射する放射輝度 $L_{i}$ は、単位投影面積当たり・単位立体角当たりの放射束で、式(2.1a)です.したがって、法線が $d\omega'$ を向いている微小面の入射束は、

d\Phi_{i} = fL_{i} \cos\psi' P(\alpha) d\omega' dA d\omega_{i}. \tag{2.6}

微小面の表面での反射は完全反射ではないので,反射する量はフレネル反射率 $F(\psi',\hat{n})$ によって決まります.$\hat{n}$ は屈折率で複素数です.

F(\psi',\hat{n}) d\Phi_{i}. \tag{2.7}

隣接する微小面は、他の微小面からの反射束,または他の微小面への入射束を遮ることがあります(マスキング,シャドウイング).係数 $G(\psi_{p},\theta_{p})$ は、実際に微小面から反射された放射束の割合です.これを幾何減衰係数といいます.したがって、領域の微小面から鏡面反射された放射束は、

\begin{align}
d\Phi_{r} &= G(\psi_{p}, \theta_{p}) F(\psi',\hat{n}) d\Phi_{i} \\[2ex]
          &= fL_{i} G(\psi_{p}, \theta_{p}) F(\psi',\hat{n}) \cos\psi' P(\alpha) d\omega' dA d\omega_{i} \tag{2.7}
\end{align}

tuto-specular-brdf-torrance-sparrow-spherical-triangle.png

反射された放射輝度 $dL_{r,s}(\psi;\theta,\phi)$ は、単位投影面当たり,単位立体角当たりの鏡面反射された放射束を表します.式(2.1b)から,単位立体角 $d\omega_{r}$ における領域 $dA$ の要素から鏡面反射する割合は $dL_{r, s}$ で表されます.

d\Phi_{r} = dL_{r,s}(\psi;\theta,\phi) \cos\theta dA d\omega_{r} \tag{2.9}

式(2.8)と(2.9)の放射束は等しく,微小面の法線方向の立体角 $d\omega'$ は以下の関係となっています.(定数 4 は Walter モデルで説明します)

d\omega' = d\omega_{r}/4 \cos\psi'. \tag{2.10}

そのため

dL_{r,s}(\psi;\theta,\phi) = (fL_{i} d\omega_{i} / 4) F(\psi',\hat{n})[G(\psi_{p},\theta_{p})/\cos\theta] P(\alpha). \tag{2.11}

式(2.3)(2.4)(2.5)を組み合わせると

dL_{r}(\psi;\theta,\phi) = \frac{fL_{i} d\omega_{i}}{4} \cdot 
F(\psi',\hat{n}) \cdot 
\frac{G(\psi_{p},\theta_{p})}{\cos\theta} \cdot 
e^{-c^2\alpha^2} + 
aL_{i} \cos\psi. \tag{2.12}

定数 $4$ はパラメータに依存しない項に,$\cos\theta$ は $G$ 項にまとめています.

この微小面による反射モデルを Blinn が簡潔にして再定義しました.BRDF $f_r$ は

f_{r} = \rho_r \frac{DGF}{\vec{n}\cdot\vec{\omega}_{o}}

ここで

  • $\rho_{r}$ は反射率です.
  • $D$ は微小面の法線分布関数です.
  • $G$ は幾何減衰項です.
  • $F$ はフレネル項です.

微小面法線分布関数には余弦分布関数,ガウス分布関数,そして Trowbridge-Reitz 分布関数があります.Trowbridge-Reitz については,Walter モデルの GGX で説明しています.ガウス分布の $\alpha$ は $\cos^{-1(\vec{n}\cdot \vec{h})}$ となっています.

幾何減衰項はV字溝によるマスキング・シャドウイングで,詳しくは「拡散BRDF」で説明しています.式で表すと

G = \min\left\{1,
  \frac{2(\vec{n}\cdot \vec{h})(\vec{n}\cdot \vec{\omega}_{o})}{(\vec{\omega}_{o}\cdot \vec{h})},
  \frac{2(\vec{n}\cdot \vec{h})(\vec{n}\cdot \vec{\omega}_{i})}{(\vec{\omega}_{i}\cdot \vec{h})}
  \right\}.

$F$ はフレネル項で,微小面がどの程度反射するかを表します.

F = \frac{1}{2} \frac{(g-c)^2}{(g+c)^2} \left\{
  1 + \frac{\left(c(g+c)-1\right)^2}{\left(c(g-c)+1\right)^2}
  \right\},

ここで:

\begin{align}
c = \cos(\theta) = \vec{\omega}_{o}\cdot \vec{h} \\[2ex]
g^2 = n^2+c^2-1
\end{align}.

法線方向(つまり $\theta = 0$ なので $c=1, g=n$)のとき,次のようになります.

F_{0} = \left\{\frac{n-1}{n+1}\right\}^2.

これから $n$ を求めると

n = \frac{1+\sqrt{F_{0}}}{1-\sqrt{F_{0}}}

視点方向と微小面法線との角度が大きくなると,視点方向から見える微小面の領域が増えるので BRDF $f_r$ の式では $\vec{n}\cdot \vec{\omega}_o$ を除算して調整しています.

微小面による反射モデルはこの Torrance-Sparrow モデルがベースになっており,微小面法線分布関数 $D$,幾何減衰項 $G$,フレネル項 $F$ を使用した反射モデルが一般的になっています.

Cook-Torrance モデル

Cook-Torrance モデルは Torrance-Sparrow モデルの法線分布項 $D$ を Beckmann 分布に変更したモデルです.また,分光反射率を用いることで,反射光の波長ごとのエネルギー分布を算出し,光源色とは異なる色として観測される金属面で鏡面反射を正確に表現できるようにしました.

表面の粗さは微小面の集合で構成されていると仮定し,$\vec{h}$ 方向の法線に向いている微小面が光源 $\vec{\omega_{i}}$ から視点 $\vec{\omega_{o}}$ に鏡面反射する放射エネルギーに寄与すると考えます.したがって,鏡面 BRDF $f_{r,s}$ は

f_{r,s} = \frac{F}{\pi} \frac{DG}{(\vec{n}\cdot \vec{\omega}_i)(\vec{n}\cdot \vec{\omega}_o)}.

$F$(フレネル項),$G$ (幾何減衰項)は Torrance-Sparrow モデルと変わりません.

微小面法線分布関数 $D$ は,$\vec{h}$ 方向に微小面がどれだけ向いているかの割合を表し,Cook-Torrance は Beckmann 分布を使いました.

D = \frac{1}{m^2\cos^4\alpha} e^{-(\tan\alpha/m)^2}

Beckmann 分布は,向きの2乗平均 $m$ を大きくすると表面が粗くなり,微小面の向きが広くなります.一方,$m$ を小さくすると,滑らかな表面となり,比較的狭い分布を生成します.その結果,鏡面反射ハイライトが鋭くなります.$\alpha$ は法線と中間ベクトルのなす角度で $\vec{n}\cdot \vec{h}$ です.

Walter モデル

Walter モデルは Cook-Torrance の微小面法線分布関数 $D$ に GGX を使いました.また,粗い表面の透過する材質をシミュレーションするために,BSDF モデルを導出しました.

BSDF(Bidirectional Scattering Distribution Function) は,双方向散乱分布関数のことで,光が表面からどのように散乱するかを表す関数です.BSDFは双方向反射率分布関数(Bidirectional Reflectance Distribution Function: BRDF)と双方向透過率分布関数(Bidirectional Transmittance Distribution Function: BTDF)の合計です.

f_{s}(\vec{\omega}_i,\vec{\omega}_o,\vec{n}) = f_{r}(\vec{\omega}_i,\vec{\omega}_o,\vec{n}) + f_{t}(\vec{\omega}_i,\vec{\omega}_o,\vec{n}).

ここで $f_{s}$ は BSDF,$f_{r}$ は BRDF,$f_{t}$ は BTDF です.今回は BTDF についての説明は割愛します.

一般的な微小面の鏡面BSDFは次のようになります.

f_{s}^m(\vec{\omega}_i,\vec{\omega}_o,\vec{m}) = \rho \frac{\delta \omega_{o}(s,\vec{\omega}_o)}{|\vec{\omega}_o\cdot\vec{m}|}.

$\delta \omega_{o}(s,\vec{\omega}_o)$ はディラックのデルタ関数で $s = \vec{\omega}_o$ のとき,値が無限大となり,それ以外の場合はゼロになる関数です.Cook-Torrance モデルと同様に微小面法線が中間ベクトルに向いている微小面が寄与していると考えます.なので,デルタ関数を使うと

f_{s}^m(\vec{\omega}_i,\vec{\omega}_o,\vec{m}) = \rho \frac{\delta \omega_{m}(\vec{h},\vec{m})}{|\vec{\omega}_o\cdot \vec{m}|}
\left\|\frac{\partial\omega_h}{\partial\omega_o}\right\|.

$\left|\frac{\partial\omega_h}{\partial\omega_o}\right|$ はヤコビアンと呼ばれるもので,変数変換を表します.これは視線方向の微小立体角 $d\omega_o$ を中間ベクトル方向の微小立体角 $d\omega_h$ に変換する係数になり,次のような極限になります.

\left\|\frac{\partial\omega_h}{\partial\omega_o}\right\| = \lim_{d\omega_{o}\to 0} \frac{d\omega_{h}}{d\omega_{o}}.

立体角は単位球面上の面積に対応していて,その領域を平面として考えると変換前と変換後の立体角の比は,変換前と変換後の微小面面積の比と同じになります.中間ベクトル $\vec{\omega} _h$ は入射ベクトル $\vec{\omega} _i$ と $\vec{\omega} _o$ を加算して求めます.

\vec{\omega}_h = \vec{\omega}_i+\vec{\omega}_o

tuto-specular-brdf-jacobian-002.png

中間ベクトルを正規化したものを $\vec{h}$ とすると

\vec{h} = \frac{\vec{\omega}_h}{\|\vec{\omega}_h\|}.

ここで立体角を面積として考えます.下図のように微小立体角 $d\omega_o$ は入射ベクトル $\vec{\omega_i}$ 移動した先の単位球上にあるので,これを元の単位球上に投影します.まずは,微小立体角 $d\omega_o$ に出射ベクトルと入射ベクトルの余弦つまり,内積の値を乗算します.

|\vec{\omega}_o\cdot \vec{\omega}_h| d\omega_o.

tuto-specular-brdf-jacobian-003.png

次に中間ベクトルの長さの二乗で投影した微小立体角 $\vec{\omega}_h$ が求まります.

d\omega_h = \frac{|\vec{\omega}_o\cdot \vec{\omega}_h| d\omega_o}{\|\vec{\omega}_h\|^2}.

tuto-specular-brdf-jacobian-004.png

ヤコビアンは変換前と変換後の比だったので

\left\|\frac{\partial\omega_h}{\partial\omega_o}\right\| = 
\frac{d\omega_h}{d\omega_o} = 
\frac{|\vec{\omega}_o\cdot \vec{\omega}_h|d\omega_o}{\|\vec{\omega}_h\|^2}\cdot\frac{1}{d\omega_o} = 
\frac{|\vec{\omega}_o\cdot \vec{\omega}_h|}{\|\vec{\omega}_h\|^2}.

中間ベクトルの長さ $|\vec{\omega_h}|$ は入射ベクトル $\vec{\omega} _i$,出射ベクトル $\vec{\omega}_o$,中間ベクトルから計算でき,内積の性質を使うと

\|\vec{\omega}_h\| = 
|\vec{\omega}_h \cdot (\vec{\omega}_i+\vec{\omega}_o)| = 
|\vec{\omega}_i\cdot\vec{\omega}_h| + |\vec{\omega}_o\cdot\vec{\omega}_h|.

$|\vec{\omega} _i\cdot\vec{\omega} _h|$ と $|\vec{\omega} _o\cdot\vec{\omega} _h|$ は同値なので

\|\vec{\omega}_h\| = 2|\vec{\omega}_i\cdot\vec{\omega}_h| = 2|\vec{\omega}_o\cdot\vec{\omega}_h|.

中間ベクトルの長さの二乗 $|\vec{\omega}_h|^2$ は

\|\vec{\omega}_h\|^2 = 
4|\vec{\omega}_i\cdot\vec{\omega}_h|^2 = 
4|\vec{\omega}_o\cdot\vec{\omega}_h|^2 = 
4|\vec{\omega}_i\cdot\vec{\omega}_h||\vec{\omega}_o\cdot\vec{\omega}_h|.

最終的にヤコビアンは

\left\|\frac{\partial\omega_h}{\partial\omega_o}\right\| = 
\frac{|\vec{\omega}_o\cdot \vec{\omega}_h|}{\|\vec{\omega}_h\|^2} = 
\frac{|\vec{\omega}_o\cdot \vec{\omega}_h|}{4|\vec{\omega}_o\cdot\vec{\omega}_h|^2} = 
\frac{1}{4|\vec{\omega}_o\cdot \vec{\omega}_h|} = \frac{1}{4|\vec{\omega}_i\cdot\vec{\omega}_h|}.

反射では $\rho$ はフレネル係数 $F$ と同じなので

f_{r}^m(\vec{\omega}_i,\vec{\omega}_o,\vec{m}) =
 F(\vec{\omega}_i,\vec{m}) \frac{\delta \omega_{m}(\vec{h},\vec{m})}{4(\vec{\omega}_i\cdot \vec{h})^2}.

$|\vec{\omega_i} \cdot \vec{h}|$ が減少すると,$f_{r}^m$ が増加し,オフスペキュラーを実現します.

表面に入射してきた放射照度と表面から出射する放射輝度は,入射角・出射角と表面の法線に基づいていますが,微小面法線が入射ベクトル方向を向いていたら通常より多く入射してくると考えられます.逆に反対の方向を向いていたら通常より少なくなるわけです.これは出射に関しても同じことがいえます.この補正をするために入射ベクトル・出射ベクトルと表面法線・微小面法線の内積の割合を乗算します.

f_{r}^m(\vec{\omega}_i,\vec{\omega}_o,\vec{m}) = 
F(\vec{\omega}_i,\vec{m}) \frac{\delta \omega_{m}(\vec{h},\vec{m})}{4(\vec{\omega}_i\cdot \vec{h})(\vec{\omega}_o\cdot \vec{h})}
\frac{(\vec{\omega}_i\cdot\vec{h})}{|\vec{\omega}_i\cdot\vec{n}|}
\frac{(\vec{\omega}_o\cdot\vec{h})}{|\vec{\omega}_o\cdot\vec{n}|}.

整理すると

f_{r}^m(\vec{\omega}_i,\vec{\omega}_o,\vec{m}) = 
F(\vec{\omega}_i,\vec{m}) 
\frac{\delta \omega_{m}(\vec{h},\vec{m})}{4|\vec{\omega}_i\cdot\vec{n}||\vec{\omega}_o\cdot \vec{n}|}.

上の式をもとに鏡面BRDFは次のようになります.

f_{r}(\vec{\omega}_i,\vec{\omega}_o,\vec{n}) = 
\frac{F(\vec{\omega}_i,\vec{h}) G(\vec{\omega}_i,\vec{\omega}_o,\vec{h}) D(\vec{h})}{4|\vec{\omega}_i\cdot \vec{n}||\vec{\omega}_o\cdot \vec{n}|}.

Cook-Torrance モデルとは分母の $\pi$ が 4 に変わっていて,それ以外の係数は同じです.$D$ は微小面法線分布関数,$G$ は幾何項,$F$ はフレネル項です.フレネル項の計算は Torrance-Sparrow モデルと変わりません.

F(\vec{\omega}_i,\vec{m}) = \frac{1}{2} \frac{(g-c)^2}{(g+c)^2}\left(
  1+\frac{(c(g+c)-1)^2}{(c(g-c)+1)^2}
  \right)\quad \text{where} \quad g = 
  \sqrt{\frac{\eta_{t}^2}{\eta_{i}^2}-1+c^2}
  \quad \text{and} \quad c = |\vec{\omega}_i\cdot \vec{m}|.

幾何項は Smith のシャドウイング・マスキングの近似式を使います.このシャドウイング・マスキング関数は双方向シャドーイングマスキング関数を2つの単方向シャドウイング項 $G_{1}$ に分離して,入射方向,出射方向と微小面法線の2つを計算しています.これを Smith Separate Masking Function といったりします.また $G_{1}$ を Smith Masking Function ともいいます.

G(\vec{\omega}_i,\vec{\omega}_o,\vec{m}) \approx G_{1}(\vec{\omega}_i,\vec{m}) G_{1}(\vec{\omega}_o,\vec{m}).

Torrance-Sparrow や Cook-Torrance のモデルではV字型溝によるシャドウイング関数でしたが,このシャドウイング関数は微小面の向きを高さと勾配で表します.具体的には微小面の高さと勾配をもとに視線方向に光線を飛ばして,どれくらい遮蔽されているかを計算します.

$G_{1}$ は微小面法線分布関数 $D$ に依存するので,この後に $D$ と一緒に明記します.任意の微小面法線分布関数 $D$ から $G_{1}$ を導出する方法が論文に書かれています.これについて簡単な説明を最後の方でします.

微小面法線分布関数 $D$ は,Beckmann, Phong, GGX の3種類が論文に記述されています.Phong と Beckmann は他のモデルでも使われていますが,BRDF の分母が異なっているので,式も他のモデルと変わっています.それぞれの分布関数と,それに対応したシャドウイング関数は次のようになっています.

$\chi^{+}(a)$ は $a > 0$ なら 1 を,そうでないなら 0 となる関数です.(ヘビサイド関数)

Beckmann 分布

D(\vec{m}) = \frac{\chi^{+}(\vec{m}\cdot \vec{n})}{\pi\alpha_{b}^2\cos^4\theta_{m}}
e^{\frac{-\tan^2\theta_{m}}{\alpha_{b}^2}}.
G_{1}(\vec{v},\vec{m}) = \chi^{+}\left(\frac{\vec{v}\cdot \vec{m}}{\vec{v}\cdot\vec{n}}\right)
\frac{2}{1+erf(a) + \frac{1}{a\sqrt{\pi}}e^{-a^2}}.

上記の式を有理数近似したものが以下の式です.

D(\vec{m}) \approx \chi^{+}\left(\frac{\vec{v}\cdot \vec{m}}{\vec{v}\cdot \vec{n}}\right)\begin{cases}
\frac{3.535\alpha+2.181\alpha^2}{1+2.276\alpha+2.577\alpha^2} && \text{if $\alpha < 1.6$} \\[2ex]
1 && \text{otherwise} \end{cases}

Phong 分布

D(\vec{m}) = \frac{\chi^{+}(\vec{m}\cdot \vec{n})} {\alpha_{p}+2}{2\pi}(\cos\theta_{m})^{\alpha_{p}}.

Beckmann 分布と同じような分布にする場合は $\alpha_{p} = \alpha_{b}^{-2}-2$ を使います.幾何項は Beckmann の幾何項の式を使って $\alpha = \sqrt{0.5\alpha_{p}+1}/(\tan\theta_{m})$ を使います.

GGX 分布

GGX 分布は Beckmann 分布や Phong 分布よりも強いテールを持っていることが特徴です.これと同じ分布関数を Trowbridge と Retiz がすでに発表していたので,Trowbrdige-Retiz 分布とも言います.

D(\vec{m}) = \frac{\alpha_{g}^2 \chi^{+}(\vec{m}\cdot \vec{n})}{
  \pi\cos^4(\theta_{m})(\alpha_{g}^2+\tan^2(\theta_{m})^2}.
G_{1}(\vec{v},\vec{m}) = \chi^{+}\left(\frac{\vec{v}\cdot \vec{m}}{\vec{v}\cdot \vec{n}}\right) \frac{2}{
  1+\sqrt{1+\alpha_{p}^2\tan(\theta_m)}}.

Phong モデル

Phong モデルは経験則に基づく反射モデルです.パラメータは物理的な意味を持っていません.ですが,計算コストが低くパラメータを調整すればそれらしい反射をしてくれるので,物理ベースレンダリング以前の鏡面反射モデルでよく使われていました.また,Phong モデルをベースに拡張された反射モデルも多く存在します.

従来の Phong モデルを BRDF として使用するとエネルギー保存則が満たされていないので,出射される放射エネルギーが入射した放射エネルギーよりも多くなることがあります.そのため,エネルギー保存則を満たすために,Phong モデルを Lewis が正規化しました.正規化されたこのモデルを正規化Phongモデルとかコサインローブモデルと言われることがあります.

f_{r,s} = \rho_s \frac{n+2}{2\pi} \cos^n \alpha.

$n$ はローブの幅を制御する係数で,$\alpha$ は入射方向に対する完全鏡面反射方向と反射方向との間の角度です.ローブの幅が広いと鈍い反射となりハイライトが大きくなります.反対に幅が狭いと鋭い反射となりハイライトが小さくなります.

Blinn が Phong モデルを改良した Blinn-Phong モデルはより現実に近い鏡面反射になります.これは Phong モデルでの完全鏡面反射方向を中間ベクトルにしたものです.この Blinn-Phong モデルも BRDF として使えるように正規化する必要があります.

f_{r,s} = \rho_s \frac{(n+2)(n+4)}{
  8\pi\left(2^{-\frac{n}{2}}+n\right)
}.

Phong モデルと Blinn-Phong モデルの違いは以下が参考になります.

正規化の計算は以下が参考になると思います.

異方性反射モデル

これまでの反射モデルは,物体表面の法線方向を軸に回転させても,反射特性が変わらないものでした.このような反射を等方性反射(isotropic reflection)といいます.通常は,入射方向・出射方向と物体表面の法線方向に対する傾きを表す角度($\theta$)と,法線方向を軸とした回転角度($\phi$)が必要です.このような反射を異方性反射(anisotropic reflection)といいます.

Ward モデル

Ward モデルは Torrance-Sparrow モデルの法線分布項 $D$ を表面の粗さ係数を対象軸に平行な方向と垂直な方向に別々に定義できるように拡張しました.これにより,法線を軸に観測面を回転させたときの見え方の変化を実現できます.

\begin{align}
f_{r,s} &= \rho_{s}
\frac{1}{\sqrt{
  (\vec{\omega}_i\cdot \vec{n})
  (\vec{\omega}_o\cdot \vec{n})
  }}
\frac{e^c}{4\pi\alpha_{x}\alpha_{y}} \\[2ex]
c &= -\tan^2 \theta_h\left(
  \frac{\cos^2\phi_h}{\alpha_{x}^2} +
  \frac{\sin^2\phi_h}{\alpha_{y}^2}
  \right).
\end{align}

$\theta_h$ は中間ベクトルと法線ベクトルの間の角度です.$\phi_h$ は中間ベクトルを表面に投影したときの方位角です.$\alpha_x$ と $\alpha_y$ は鏡面ローブ幅の係数で,$\alpha_x = \alpha_y$ のとき,等方性反射になります.

上記の近似式は次のようになっています.三角関数の部分はベクトルの演算に置き換わっています.

c = -2
\left\{
  \left(\frac{\vec{h}\cdot\vec{x}}{\alpha_{x}}\right)^2+
  \left(\frac{\vec{h}\cdot\vec{y}}{\alpha_{y}}\right)^2
\right\}
\frac{1}{1+\vec{h}\cdot\vec{n}}.

$\vec{x}$ と $\vec{y}$ は法線ベクトルに対する接ベクトル,従法線ベクトルです.Walter は上記の近似式よりも正確な近似式を提示しました.

c = -\left\{
  \left(\frac{\vec{h}\cdot\vec{x}}{\alpha_{x}}\right)^2+
  \left(\frac{\vec{h}\cdot\vec{y}}{\alpha_{y}}\right)^2
\right\}
\frac{1}{(\vec{h}\cdot\vec{n})^2}

Ashikhmin モデル

Phong モデルを異方性反射に対応した反射モデルです.具体的にはエネルギー保存則と相反性に従っていること,異方性反射であること,パラメータが直感的であること,フレネル効果を考慮しているといった特徴があります.

\begin{align}
f_{r,s} &= \rho_s \frac{\sqrt{(\alpha_x+1)(\alpha_y+1)}
}{
  8\pi
}
\frac{
  (\vec{n}\cdot \vec{h})^c
}{
  (\vec{h}\cdot \vec{\omega})\max((\vec{n}\cdot \vec{\omega}_i),(\vec{n}\cdot \vec{\omega}_o))
}
F(\vec{\omega}\cdot \vec{h}) \\[2ex]
c &= \frac{
  (\alpha_x(\vec{x}\cdot\vec{h})^2+\alpha_y(\vec{y}\cdot\vec{h})^2)
}{
  (1-(\vec{h}\cdot\vec{n})^2)
}.
\end{align}

$\alpha_x$ と $\alpha_y$ は Ward モデルと同様に鏡面ローブ幅の係数で,$\alpha_x = \alpha_y$ のとき,等方性反射になります.$\vec{\omega}$ は $\vec{\omega}_i$ または $\vec{\omega}_o$ のどちらでも使えます.$F$ はフレネル項で,Schlick の近似式を使います.

F(\vec{\omega}\cdot \vec{h}) = \rho_s + (1-\rho_s)(1-(\vec{\omega}\cdot \vec{h}))^5.

余談になりますが,Ashikhmin の拡散反射モデルはランバート反射にフレネル効果を追加したものになっています.

f_{r,d} = \frac{\rho_d}{\pi} \frac{28}{23}(1-\rho_s)\left(
  1-\left(1-\frac{(\vec{\omega}_i\cdot\vec{n})}{2}\right)^5
\right)\left(
  1-\left(1-\frac{(\vec{\omega}_o\cdot\vec{n})}{2}\right)^5
\right).

$\frac{28}{23}$ はエネルギー保存則を満たすための正規化係数です.

BRDF は拡散BRDFと鏡面反射BRDFを加算したものです.

f_r = f_{r,d} + f_{r,s}.

Schlick モデル

物理ベースモデルでは,光の正確な物理現象を再現することができます.しかし,計算も複雑になってしまいます.そこで,Schlick モデルでは正しい計算と同じ結果になるような近似式を使って,微小面モデルの計算を軽くし,異方性反射に対応したモデルです.このモデルでは以下の3つのパラメータ(反射係数,表面の粗さ,等方性係数)を使います.

  • $F_0 \in [0,1]$ : 表面に垂直入射したときの反射係数.
  • $\sigma \in [0,1]$ : 粗さの要素.$\sigma = 0$ で完全に滑らかな鏡面で,$\sigma = 1$ で非常に粗いランバート面になります.
  • $\psi \in [0,1]$ : 等方性係数.$\psi = 0$ で完全な異方性,$\phi = 1$ で等方性になります.
f_r = S(u)\left\{
  \frac{d}{\pi}+
  gD(t,v_o,v_i,w) +
  s f_{r,s}(x,\vec{\omega}_o,\vec{\omega}_i)
\right\}.

Schlick モデルの BRDF は拡散反射と光沢反射,鏡面反射の和となっています.それぞれ $d, g, s$ が係数となっています.その他の各係数は

\begin{align}
u &= \vec{\omega}_o\cdot \vec{h} \\[2ex]
t &= \vec{h}\cdot \vec{n} \\[2ex]
v_o &= \vec{\omega}_o \cdot \vec{n} \\[2ex]
v_i &= \vec{\omega}_i \cdot \vec{n} \\[2ex]
w &= \vec{T}\cdot \frac{
  \vec{h} - (\vec{h}\cdot\vec{n})\vec{n}
}{
  |\vec{h} - (\vec{h}\cdot\vec{n})\vec{n}|
}
\end{align}

$\vec{T}$ は異方性物質における傷(scratch)方向を向いた接ベクトルです.$S(u)$ は鏡面反射成分で,Schlick のフレネル近似で求まります.

S(u) = F_0 + (1-F_0)(1-u)^5.

$f_{r,s}$ は次のような完全鏡面反射です.

f_{r,s} = 2 \rho_s
 \delta(\sin^2(\theta_o)-\sin^2(\theta_i))
 \delta(\phi_i-\phi_o\pm\pi). 

$D$ 項は方向成分で,光沢反射を表します.一般的な鏡面反射はこの光沢反射と完全鏡面反射を足したものと考えることができます.

D(t,v_o,v_i,w) = \frac{1}{4\pi v_o v_i} Z(t) A(w).

ここで

Z(t) = \frac{\sigma}{(1+\sigma t^2-t^2)^2} \quad and \quad 
A(w) = \sqrt{\frac{\psi}{\psi^2-\psi^2w^2+w^2}}.

次に,D項に幾何減衰係数を追加してシャドウイング・マスキングの効果を適用します.

D(t,v_o,v_i,w) = \frac{G(v_o)G(v_i)}{4\pi v_o v_i} Z(t) A(w) +
\frac{1-G(v_o)G(v_i)}{4\pi v_o v_i}.

$G$ は Smith の近似式を使って

G(v_o) = \frac{v_o}{
  \sigma - \sigma v_o + v_o
} \quad \text{and} \quad 
G(v_i) = \frac{v_i}{
  \sigma - \sigma v_i + v_i
}.

拡散反射係数 $d$,光沢反射係数 $g$,完全鏡面反射係数 $s$ は表面の粗さ $\sigma$ に基づいて次のように決定することが提案されています.

\begin{align}
g &= 4\sigma(1-\sigma) \\[2ex]
d &= \begin{cases} 0 && \text{if $\sigma < 0.5$} \\[2ex]
 1-g && \text{otherwise} \end{cases} \\[2ex]
s &= \begin{cases} 1-g && \text{if $\sigma < 0.5$} \\[2ex]
0 && \text{otherwise} \end{cases} \\[2ex]
\end{align}

上記の式に従えば $g+d+s=1$ となり,エネルギー保存則を満たします.

Disney モデル

Walter の微小面モデルをもとに改良したモデルで,次のような原則に従っています.

  • パラメータは直感的である
  • 可能な限り少ないパラメータにする
  • パラメータは [0,1] の範囲にする
  • 意図的にパラメータを範囲外の値にしても,物理的に問題がなければ許容する
  • パラメータのすべての組み合わせにおいて,ロバストであること

このため,完全に物理ベースというわけではないのですが,ほぼ物理ベースといえるでしょう.

Disney モデルの BRDF は他のモデルと同様に,拡散 BRDF と鏡面 BRDF の和です.Disney モデルの拡散BRDFについては「拡散BRDF」で説明しています.

パラメータは次のようになっています.

  • baseColor: 表面の色です.
  • subsurface: 表面下の近似を使って拡散反射を制御します.
  • metallic: 金属度合いです.0 は誘電体,1 が金属で,中間はそれらの線形合成になります.
  • specular: 入射鏡面反射で,屈折率の代わりです.
  • specularTint: 鏡面反射を表面の色に近づける度合いです.
  • roughness: 表面の粗さです.
  • anisotropic: 異方性の度合いです.0 のとき,等方性になります.
  • sheen: 布用のグレージング角に対する追加成分です.
  • sheenTint: shenn を表面の色に近づける度合いです.
  • clearcoat: 2層目のクリアコートの強さです.
  • clearcoatGloss: クリアコートの光沢度です.0 がサテン風,1 がグロス風になります.

tuto-specular-brdf-disney-parameters.png

Brent Burley, "Physically-Based Rendering at Disney", p. 13 より引用

法線分布関数 $D$ は GTR を使用します.これは Generalized-Trowbridge-Reitz の略で,GGX (Trowbridge-Reitz)分布をより一般的にしたものです.

D_{GTR} = \frac{c}{(\alpha^2\cos^2\theta_h + \sin^2\theta_h)^{\gamma}}.

$c$ はスケール定数です.$\gamma = 2$ とすると GGX 分布と同じになります.

tuto-specular-brdf-disney-gtr-distribution.png

Brent Burley, "Physically-Based Rendering at Disney", p. 15 より引用

Disney モデルの反射は2層となっていて,2層目はクリアコートを表しています.1層目は $\gamma = 2$ の GGX 分布を使用し,2層目は $\gamma = 1$ を使っていて,これは等方性で非金属となります.

D_{GTR_{1}}(\theta_h) = \frac{
  \alpha^2-1
}{
  \pi \log \alpha^2
} \frac{1}{
  (1+(\alpha^2-1)\cos^2\theta_h)
}.
D_{GTR_{2aniso}} = \frac{1}{\pi}
  \frac{1}{\alpha_x \alpha_y}
  \frac{1}{
    \left(
      (\vec{h}\cdot\vec{x})^2 / \alpha_x^2 +
      (\vec{h}\cdot\vec{y})^2 / \alpha_y^2 +
      (\vec{h}\cdot\vec{n})^2
    \right)^2
  }.

表面の粗さを表す $\alpha$ は,roughness を2乗した値となっています.これは,粗さによる知覚的な変化がより線形になるからです.

フレネル項 $F$ は Schlick の近似式を使います.

F_{schlick} = F_0 + (1-F_0)(1-\cos\theta_d)^5.

幾何減衰項 $G$ は Walter モデルの $G$ を使っていますが,表面の粗さ $\alpha_g$ は,$(0.5 + \text{roughness}/2)^2$ としています.

G_{GGX}(\vec{v}) = \frac{1}{
  \left(
    (\vec{v}\cdot \vec{n}) +
    \sqrt{
      \alpha_g^2 + 
      (\vec{v}\cdot \vec{n})^2 - 
      \alpha_g^2 (\vec{v}\cdot \vec{n})^2
    }
  \right)
}.
G_{GGX_{aniso}} = \frac{1}{
  \left(
    (\vec{v}\cdot \vec{n}) +
    \sqrt{
      \left(
        (\vec{v}\cdot \vec{x}) a_x
      \right)^2 +
      \left(
        (\vec{v}\cdot \vec{y}) a_y
      \right)^2 +
      (\vec{v}\cdot \vec{n})^2
    }
  \right)
}.

Disney モデルの BRDF $f_r$ は拡散 BRDF,鏡面 BRDF,クリアコート BRDF の和になります.

f_{r} = f_{r,d} + f_{r,s} + f_{r,clearcoat}

鏡面 BRDF $f_r$ は次のようになります:

f_{r,s} = G_s F_s D_s.

各係数は:

\begin{align}
\text{mix}(a,b,t) &= a + (b-a)t \\[2ex]
\text{luminance} &= 0.3\cdot \text{baseColor}_r + 0.6\cdot \text{baseColor}_g + 0.1\cdot \text{baseColor}_b \\[2ex]
\text{Ctint} &= \frac{\text{baseColor}_{rgb}}{\text{luminance}} \\[2ex]
\text{Cspec} &= \text{mix}(\text{specular}\cdot 0.08 \cdot \text{mix}(1, \text{Ctint}, \text{specularTint}) \\[2ex]
\text{aspect} &= \sqrt{1-\text{anisotropic}\cdot 0.9} \\[2ex]
a_x &= \max(0.001, \frac{\text{roughness}^2}{\text{aspect}}) \\[2ex]
a_y &= \max(0.001, \text{roughness}^2 \cdot \text{aspect}) \\[2ex]
D_s &= D_{GTR_{2aniso}}(
  (\vec{h}\cdot\vec{n}),
  (\vec{x}\cdot\vec{h}),
  (\vec{y}\cdot\vec{h}),
  a_x, a_y) \\[2ex]
G_s &= G_{GGX_{aniso}}(
  (\vec{\omega}_i,\vec{n}),
  (\vec{\omega}_i,\vec{x}),
  (\vec{\omega}_i,\vec{y}),
  a_x, a_y) \cdot
  G_{GGX_{aniso}}(
    (\vec{\omega}_o,\vec{n}),
    (\vec{\omega}_o,\vec{x}),
    (\vec{\omega}_o,\vec{y}),
    a_x, a_y) \\[2ex]
F_h &= F_{schlick}(\vec{\omega}_i,\vec{h}) \\[2ex]
F_{s} &= \text{mix}(\text{Cspec}, 1, F_h) \\[2ex]
\end{align}

拡散 BRDF $f_d$ は次のようになります:

f_{r,d} = \left(
  \frac{1}{\pi} \cdot
  \text{mix}(F_d, \text{ss}, \text{subsurface}) \cdot
  \text{baseColor}_{rgb} + CF_{sheen}
\right) \cdot (1-\text{metallic}).

各係数は:

\begin{align}
\text{Csheen} &= \text{mix}(1, \text{Ctint}, \text{sheenTint}) \\[2ex]
F_i &= F_{schlick}(\vec{\omega}_i\cdot\vec{n}) \\[2ex]
F_o &= F_{schlick}(\vec{\omega}_o\cdot\vec{n}) \\[2ex]
F_{d90} &= 0.5 + 2\cdot (\vec{\omega}_i\cdot\vec{h})(\vec{\omega}_o\cdot\vec{h})\cdot \text{roughness} \\[2ex]
F_{d} &= \text{mix}(1,F_{d90},F_i)\cdot \text{mix}(1,F_{d90},F_o) \\[2ex]
F_{s90} &= (\vec{\omega}_i\cdot\vec{h})(\vec{\omega}_o\cdot\vec{h})\cdot \text{roughness} \\[2ex]
F_{ss} &= \text{mix}(1,F_{s90},F_i)\cdot\text{mix}(1,F_{s90},F_o) \\[2ex]
\text{ss} &= 1.25\cdot \left(
  F_{ss} \left( \frac{1}{
    (\vec{\omega}_i\cdot\vec{n})(\vec{\omega}_o\cdot\vec{n})
    } - 0.5) + 0.5
  \right)
\right) \\[2ex]
CF_{sheen} &= F_h \cdot \text{Csheen} \cdot \text{sheen} \\[2ex]
\end{align}

2層目のクリアコート BRDF $f_{r,clearcoat}$ は次のようになります:

f_{r,clearcoat} = 0.25\cdot \text{clearcoat} \cdot G_r F_r D_r.

各係数は:

\begin{align}
D_r &= D_{GTR_{1}}(
  (\vec{h}\cdot\vec{n}),
  \text{mix}(0.1, 0.001, \text{clearcoatGloss})) \\[2ex]
F_r &= \text{mix}(0.04, 1, F_h) \\[2ex]
G_r &= G_{GGX}((\vec{\omega}_i\cdot\vec{n}), 0.25) \cdot 
G_{GGX}((\vec{\omega}_o\cdot\vec{n}), 0.25) \\[2ex]
\end{align}

微小面の法線分布関数

説明していない微小面の法線分布関数に ABC 分布と SGD 分布があります.詳細は[Hoffman13]や[Hanecci13]を参照してください.

Smith Shadowing-Masking Function

Separable Masking and Shadowing

単方向シャドウイング関数 $G_{1}$ を使って,入射と出射方向の双方向によるシャドウイングを計算し,それらを乗算することでシャドウイング・マスキングを計算するモデルです.定義を見ると

G(\vec{\omega}_i,\vec{\omega}_o) \approx G_1(\vec{\omega}_i) G_1(\vec{\omega}_o) = \frac{1}{1 + \Lambda(\vec{\omega}_i)}\cdot \frac{1}{1+\Lambda(\vec{\omega}_o)}.

ここで $\Lambda$ は微小面法線分布関数 $D$ によって変わる関数です.任意の $D$ に対応した $G_{1}$ の導出についてわかる範囲で説明します.

基本的な考え方ですが,まず,各微小面の開始位置における高さと勾配をそれぞれの確率分布から決定します.そうすると表面は次のような微小面で構成されています.

tuto-specular-brdf-masking-function-002.png

微小面の最初の位置における高さ $\xi$ から視線方向に光線を飛ばして隠蔽密度を調べます.

tuto-specular-brdf-masking-function-003.png

微小面の高さは2つの確率分布で表されると仮定します.$P_{1}$ は高さ $\xi$ の確率分布関数で,$P_{22}(p,q)$ は入射面に対する垂直方向と水平方向の勾配 $p$ と $q$ の確率分布関数です.

$P_{1}$ は任意の確率分布関数を指定することができ,$P_{22}$ は微小面法線分布関数 $D$ から次のように計算することができます.

P_{22}(p,q) = D(\vec{m})\cos^4\theta_m.

ここで,等方性の反射モデルと仮定すると,水平方向の勾配 $p$ は一定となるので,垂直方向の勾配 $q$ の累積分布関数 $P_{2}$ は次のように計算できます.

P_{2}(q) = \int_{-\infty}^{\infty} P_{22}(p,q)dp.

$S(\xi_0,\mu)$ は $\xi_0$ の高さを持つ微小面上の無作為に選んだ点が視線方向 $\vec{v}$ から見える確率です.次に光線に必要なパラメータを定義します.光線の開始位置は微小面上の高さ $\xi_0$ の点です.光線の勾配 $\mu$ は表面法線 $\vec{n}$ と視線方向 $\vec{v}$ から次のように求められます.

\mu = |\cot\theta_v|.

tuto-specular-brdf-masking-function-001.png

光線の開始位置から視線方向に進んだ分,表面上に投影した距離を $\tau$ とすると,$\tau$ における光線の高さは $\xi_0 + \mu\tau$ となります.区間[$\tau,\tau+\Delta\tau$]で,遮蔽された割合を $g(\tau)\Delta\tau$ とすると,$S(\xi_0,\mu)$ は

S(\xi_0,\mu) = e^{-\int_{0}^{\infty} g(\tau)d\tau}.

$\tau$ における高さと勾配がそれぞれ $\xi,q$ とすると,$\xi_0+\mu\tau > \xi$ のとき,光線は $\tau$ での高さよりも上を通過します.また,$(q-\mu)\Delta\tau > (\xi_0 + \mu\tau) - \xi$ なら,光線は区間[$\tau,\tau+\Delta\tau$]のどこかで微小面と交差するか,遮蔽されます.微小面の高さと勾配の分布がそれぞれ独立していると仮定して,$g(\tau)$ を次のように近似します.

\begin{align}
g(\tau) &= \frac{
  \int_{\mu}^{\infty} (q-\mu) P_{1}(\xi_0+\mu\tau) P_2(q) dq
}{
  \int_{-\infty}^{\xi_0+\mu\tau} P_1(\xi) d\xi
} \\[2ex]
 &= \Lambda(\mu) \frac{
   \mu P_1(\xi_0+\mu\tau)
 }{
   f(\xi_0+\mu\tau)
 }
\end{align}

$f(z)$ は高さの累積密度関数で,高さ $z$ が表面より上である確率を表します.

f(z) = \int_{-\infty}^{z} P_1(\xi) d\xi.

$\Lambda(\mu)$ は勾配の平均です.

\Lambda(\mu) = \frac{1}{\mu} \int_{\mu}^{\infty} (q-\mu) P_2(q) dq.

この $g(\tau)$ の式を使って,$S(\xi_0,\mu)$ の式を解くと

S(\xi_0,\mu) = e^{\Lambda(\mu) \log f(\xi_0)} = f(\xi_0)^{\Lambda(\mu)}.

$\xi_0$ の高さから全体にわたって積分して,可視密度の平均 $S(\mu)$ を求めると次のようになります.

S(\mu) = \int_{-\infty}^{\infty} S(\xi_0,\mu) P_1(\xi_0) d\xi_0 =
\frac{1}{1 + \Lambda}.

問題点

この計算モデルにはいくつか問題があります.

微小面のすべての高さについて積分しているので,同じ高さの微小面が複数存在する場合に,どの微小面が計算されるかで計算結果が異なり,計算誤差となります.

tuto-specular-brdf-masking-function-005.png

微小面の高さと勾配が隣の微小面に依存していないので,連続になっていません.

tuto-specular-brdf-masking-function-004.png

Height-Correlated Masking and Shadowing

微小面の高さに相関したシャドウイング・マスキング関数です.Smith Joint Masking-Shadowing Function ということもあります.これは,微小面の位置が高くなるほど,入射および出射方向に対して見える可能性が大きくなることを考慮した関数で,次のように定義されています.

G(\vec{\omega}_i,\vec{\omega}_v) = \frac{1}{1+\Lambda(\vec{\omega}_l)+\Lambda(\vec{\omega}_v)}.

この式は,入射方向と出射方向が近くなるほどシャドウイングが大きくなります.

最後に

いくつかの鏡面BRDF,そして反射モデルを紹介しました.まだまだ多くの反射モデルがあって,Lafortune モデル,He-Torrance モデル,Kajiya-Kay モデル,Marschner モデルというのがあります.他にも改良されたモデルや新しいモデルがあると思います.今回の内容について,私の知識不足,理解不足で間違っているところがあると思いますので,見つけたらご連絡していただけると助かります.少しでも参考になれば幸いです.

関連記事

参考資料