Microfacet入門（２） - ushiostarfish’s diary

スロープ空間

先のMicrofacet入門（１）でおざなりにしていた、 $P_{ 2 }$ とは何かについて、見ていきます。それを整理するために、まずはスロープ空間について考えます。

スロープ空間とは、マイクロファセット法線を変数変換してできる空間です。マイクロファセット法線を、 $\theta, \phi$ を使って表すと、

$\displaystyle { \omega _{ m }=\left( \cos { \phi } \sin { \theta } ,\sin { \theta } \sin { \phi } ,\cos { \theta } \right) \\ \omega _{ m }=\left( x_{ m },y_{ m },z_{ m } \right) }$

となり、スロープ空間は二次元であり、

$\displaystyle { \tilde { m } \left( \omega _{ m } \right) =\left( x_{ \tilde { m } },y_{ \tilde { m } } \right) }$

と表すことにすると、変数変換は、

$\displaystyle { \tilde { m } \left( \omega _{ m } \right) =\left( -\frac { x_{ m } }{ z_{ m } } ,-\frac { y_{ m } }{ z_{ m } } \right) \\ \tilde { m } \left( \omega _{ m } \right) =\left( -\frac { \cos { \phi } \sin { \theta } }{ \cos { \theta } } ,-\frac { \sin { \theta } \sin { \phi } }{ \cos { \theta } } \right) }$

と表すことができます。一見直観的でないように思われますが、 $\phi =-\pi$ である状況を考えてみましょう。

$\displaystyle { x_{ \tilde { m } }=-\frac { \cos { \left( -\pi \right) } \sin { \theta } }{ \cos { \theta } } =\frac { \sin { \theta } }{ \cos { \theta } } =\tan { \theta } \\ y_{ \tilde { m } }=-\frac { \sin { \theta } \sin { \left( -\pi \right) } }{ \cos { \theta } } =0 }$

なんと、x が $\tan \theta$ になりました。これが"スロープ空間"のゆえんであります。これには幾何学的解釈もでき、以下のようにスロープ空間がxy平面だとすると、z=1の面に投射した、と考えることができます。

f:id:ushiostarfish:20180401061608p:plain ※75 ページ Earl Hammon, Jr, GDC, "PBR Diffuse Lighting for GGX+Smith Microsurfaces" より

なぜかといえば、法線のxyは、zで割り算されているためで、二次元的には以下のような状態だからです。

f:id:ushiostarfish:20180401064853p:plain

ちなみに逆変換も可能で、

$\displaystyle { x_{ \tilde { m } }=-\frac { x_{ m } }{ z_{ m } } \\ y_{ \tilde { m } }=-\frac { y_{ m } }{ z_{ m } } \\ -x_{ \tilde { m } }z_{ m }=x_{ m }\\ -y_{ \tilde { m } }z_{ m }=y_{ m } }$

であるため、法線 $\omega _m$ は、

$\displaystyle { \omega _{ m }=\left( -x_{ \tilde { m } }z_{ m },-y_{ \tilde { m } }z_{ m },z_{ m } \right) }$

ただ、 $z_{ m }$ は未知であり、

$\displaystyle { \left| \omega _{ m } \right| =\sqrt { { x_{ m } }^{ 2 }+{ y_{ m } }^{ 2 }+{ z_{ m } }^{ 2 } } =1\\ \sqrt { { \left( x_{ \tilde { m } }z_{ m } \right) }^{ 2 }+{ \left( y_{ \tilde { m } }z_{ m } \right) }^{ 2 }+{ z_{ m } }^{ 2 } } =1\\ { \left( x_{ \tilde { m } }z_{ m } \right) }^{ 2 }+{ \left( y_{ \tilde { m } }z_{ m } \right) }^{ 2 }+{ z_{ m } }^{ 2 }=1\\ { \left( x_{ \tilde { m } } \right) }^{ 2 }+{ \left( y_{ \tilde { m } } \right) }^{ 2 }+1=\frac { 1 }{ { z_{ m } }^{ 2 } } \\ { z_{ m } }^{ 2 }=\frac { 1 }{ { \left( x_{ \tilde { m } } \right) }^{ 2 }+{ \left( y_{ \tilde { m } } \right) }^{ 2 }+1 } \\ { z_{ m } }=\frac { 1 }{ \sqrt { { \left( x_{ \tilde { m } } \right) }^{ 2 }+{ \left( y_{ \tilde { m } } \right) }^{ 2 }+1 } } }$

とすると、求めることができるため、最終的には、

$\displaystyle { \omega _{ m }=\frac { 1 }{ \sqrt { { \left( x_{ \tilde { m } } \right) }^{ 2 }+{ \left( y_{ \tilde { m } } \right) }^{ 2 }+1 } } \left( -x_{ \tilde { m } },-y_{ \tilde { m } },1 \right) }$

となります。法線が逆側の場合はどうか？というのもありますが、実際ミクロ表面は裏返らないため、問題ありません。

スロープ空間とヤコビアン

変換方法は明らかになりましたが、元の"法線分布"を引き継ぐため、この変数変換のヤコビアンを求めておきたいです。もとの変数は $\theta, \phi$ であり、ここから $x_{ \tilde { m } },y_{ \tilde { m } }$ への変換は、ヤコビ行列は、

$\displaystyle { J=\begin{pmatrix} \frac { \partial x_{ \tilde { m } } }{ \partial \theta } & \frac { \partial x_{ \tilde { m } } }{ \partial \phi } \\ \frac { \partial y_{ \tilde { m } } }{ \partial \theta } & \frac { \partial y_{ \tilde { m } } }{ \partial \phi } \end{pmatrix}=\begin{pmatrix} -\cos \left( \phi \right) \sec ^{ 2 } \left( \theta \right) & \sin \left( \phi \right) \tan \left( \theta \right) \\ -\sin \left( \phi \right) \sec ^{ 2 } \left( \theta \right) & -\tan \left( \theta \right) \cos \left( \phi \right) \end{pmatrix} }$

と、これは機械的に求め、ここからヤコビアンは、

$\displaystyle { \left| J \right| =\cos \left( \phi \right) \sec ^{ 2 } \left( \theta \right) \tan \left( \theta \right) \cos \left( \phi \right) +\sin \left( \phi \right) \tan \left( \theta \right) \sin \left( \phi \right) \sec ^{ 2 } \left( \theta \right) \\ =\cos ^{ 2 } \left( \phi \right) \sec ^{ 2 } \left( \theta \right) \tan \left( \theta \right) +\sin ^{ 2 } \left( \phi \right) \sec ^{ 2 } \left( \theta \right) \tan \left( \theta \right) \\ =\sec ^{ 2 } \left( \theta \right) \tan \left( \theta \right) \left( \cos ^{ 2 } \left( \phi \right) +\sin ^{ 2 } \left( \phi \right) \right) \\ =\frac { 1 }{ \cos ^{ 2 }{ \theta } } \frac { \sin { \theta } }{ \cos { \theta } } \\ =\frac { \sin { \theta } }{ \cos ^{ 3 }{ \theta } } }$

と求めることができます。ヤコビアンは微小な成分の拡大率であるため、以下のように積分の辻褄があったまま変数変換することができます。

$\displaystyle { \int _{ -\infty }^{ +\infty }{ \int _{ -\infty }^{ +\infty }{ f\left( x_{ \tilde { m } },y_{ \tilde { m } } \right) dx_{ \tilde { m } } } dy_{ \tilde { m } } } =\int _{ 0 }^{ 2\pi }{ \int _{ 0 }^{ \frac { \pi }{ 2 } }{ f\left( \theta ,\phi \right) \frac { \sin { \theta } }{ \cos ^{ 3 }{ \theta } } d\theta d\phi } } }$

これを使って、

$\displaystyle { 1=\int _{ \Omega }^{ }{ D\left( \omega _{ m } \right) \left| \omega _{ m }\cdot \omega _{ g } \right| d\omega _{ m } } \\ 1=\int _{ -\infty }^{ +\infty }{ \int _{ -\infty }^{ +\infty }{ { P }^{ 22 }\left( x_{ \tilde { m } },y_{ \tilde { m } } \right) dx_{ \tilde { m } } } dy_{ \tilde { m } } } }$

を満たすような ${ P }^{ 22 }$ を考えます。これには、

$\displaystyle { \int _{ \Omega }^{ }{ D\left( \omega _{ m } \right) \left| \omega _{ m }\cdot \omega _{ g } \right| d\omega _{ m } } ={ \int _{ 0 }^{ 2\pi }{ \int _{ 0 }^{ \frac { \pi }{ 2 } }{ D\left( \omega _{ m } \right) \cos { \theta } \sin { \theta } d\theta d\phi } } } }$

であり、

$\displaystyle { \int _{ -\infty }^{ +\infty }{ \int _{ -\infty }^{ +\infty }{ { P }^{ 22 }\left( x_{ \tilde { m } },y_{ \tilde { m } } \right) dx_{ \tilde { m } } } dy_{ \tilde { m } } } =\int _{ 0 }^{ 2\pi }{ \int _{ 0 }^{ \frac { \pi }{ 2 } }{ { P }^{ 22 }\left( x_{ \tilde { m } },y_{ \tilde { m } } \right) \frac { \sin { \theta } }{ \cos ^{ 3 }{ \theta } } d\theta d\phi } } }$

であることから、

$\displaystyle { { P }^{ 22 }\left( x_{ \tilde { m } },y_{ \tilde { m } } \right) \frac { \sin { \theta } }{ \cos ^{ 3 }{ \theta } } =D\left( \omega _{ m } \right) \cos { \theta } \sin { \theta } \\ { P }^{ 22 }\left( x_{ \tilde { m } },y_{ \tilde { m } } \right) =D\left( \omega _{ m } \right) \cos ^{ 4 }{ \theta } }$

のようにDから機械的に変換することができることがわかります。

具体的なP

そろそろ具体的なDを考えて、そこからPを出しておきましょう。有名な法線分布として今回は、Beckmann分布と、GGX(Trowbridge-Reitz)分布を取り上げます。まずBeckmann分布ですが、これの具体的な形は $\alpha$ をパラメーターとして、

$\displaystyle { D_{ beckmann }\left( \omega _{ m } \right) =\frac { { \chi }^{ + }\left( \omega _{ m }\cdot \omega _{ g } \right) }{ \pi { \alpha }^{ 2 }\cos ^{ 4 }{ \theta } } exp\left( -\frac { \tan ^{ 2 }{ \theta } }{ { \alpha }^{ 2 } } \right) }$

ここから ${ P }^{ 22 }$ は、

$\displaystyle { { P }^{ 22 }\left( x_{ \tilde { m } },y_{ \tilde { m } } \right) =\frac { 1 }{ \pi { \alpha }^{ 2 } } exp\left( -\frac { { x_{ \tilde { m } } }^{ 2 }+{ y_{ \tilde { m } } }^{ 2 } }{ { \alpha }^{ 2 } } \right) }$

となります。 ${ x_{ \tilde { m } } }^{ 2 }+{ y_{ \tilde { m } } }^{ 2 }$ はスロープ空間上でのゼロ地点からの距離の二乗であり、それは結局 $\tan^2 \theta$ であることを利用しています。

そう、なんとBeckmann分布はスロープ空間において、標準的な2次元ガウス関数だったのです。

続いて、GGX(Trowbridge-Reitz)分布です。

$\displaystyle { D_{ ggx }\left( \omega _{ m } \right) =\frac { { \chi }^{ + }\left( \omega _{ m }\cdot \omega _{ g } \right) }{ \pi { \alpha }^{ 2 }\cos ^{ 4 }{ \theta _{ m } } { \left( 1+\frac { \tan ^{ 2 }{ \theta _{ m } } }{ { \alpha }^{ 2 } } \right) }^{ 2 } } }$

こちらも同様に、

$\displaystyle { { P }^{ 22 }\left( x_{ \tilde { m } },y_{ \tilde { m } } \right) =\frac { 1 }{ \pi { \alpha }^{ 2 }{ \left( 1+\frac { { x_{ \tilde { m } } }^{ 2 }+{ y_{ \tilde { m } } }^{ 2 } }{ { \alpha }^{ 2 } } \right) }^{ 2 } } }$

となります。

P2

では Microfacet入門（１）での $P_{ 2 }$ は結局何かと言えば、片方を単純に積分して一次元に落としたものです。

$\displaystyle { { P }^{ 2- }\left( x_{ \tilde { m } } \right) =\int _{ -\infty }^{ \infty }{ { P }^{ 22 }\left( x_{ \tilde { m } },y_{ \tilde { m } } \right) dy_{ \tilde { m } } } }$

これは実際のところ、Beckmann分布とGGX(Trowbridge-Reitz)分布で解析的に解くことができ、

Beckmann分布の場合は、

$\displaystyle { \\ { P }^{ 2- }\left( x_{ \tilde { m } } \right) =\frac { 1 }{ \pi { \alpha }^{ 2 } } \int _{ -\infty }^{ \infty }{ exp\left( -\frac { { x_{ \tilde { m } } }^{ 2 } }{ { \alpha }^{ 2 } } -\frac { { y_{ \tilde { m } } }^{ 2 } }{ { \alpha }^{ 2 } } \right) dy_{ \tilde { m } } } \\ =\frac { 1 }{ \pi { \alpha }^{ 2 } } \int _{ -\infty }^{ \infty }{ exp\left( -\frac { { x_{ \tilde { m } } }^{ 2 } }{ { \alpha }^{ 2 } } -\frac { { y_{ \tilde { m } } }^{ 2 } }{ { \alpha }^{ 2 } } \right) dy_{ \tilde { m } } } \\ =\frac { 1 }{ \pi { \alpha }^{ 2 } } \int _{ -\infty }^{ \infty }{ exp\left( -\frac { { x_{ \tilde { m } } }^{ 2 } }{ { \alpha }^{ 2 } } -\frac { { y_{ \tilde { m } } }^{ 2 } }{ { \alpha }^{ 2 } } \right) dy_{ \tilde { m } } } \\ =\frac { 1 }{ \pi { \alpha }^{ 2 } } \int _{ -\infty }^{ \infty }{ exp\left( -\frac { { x_{ \tilde { m } } }^{ 2 } }{ { \alpha }^{ 2 } } \right) exp\left( -\frac { { y_{ \tilde { m } } }^{ 2 } }{ { \alpha }^{ 2 } } \right) dy_{ \tilde { m } } } \\ =\frac { 1 }{ \pi { \alpha }^{ 2 } } exp\left( -\frac { { x_{ \tilde { m } } }^{ 2 } }{ { \alpha }^{ 2 } } \right) \int _{ -\infty }^{ \infty }{ exp\left( -\frac { { y_{ \tilde { m } } }^{ 2 } }{ { \alpha }^{ 2 } } \right) dy_{ \tilde { m } } } \\ =\frac { 1 }{ \pi { \alpha }^{ 2 } } exp\left( -\frac { { x_{ \tilde { m } } }^{ 2 } }{ { \alpha }^{ 2 } } \right) \sqrt { \pi { \alpha }^{ 2 } } \\ =\frac { \sqrt { \pi { \alpha }^{ 2 } } }{ \pi { \alpha }^{ 2 } } exp\left( -\frac { { x_{ \tilde { m } } }^{ 2 } }{ { \alpha }^{ 2 } } \right) \\ =\frac { \sqrt { \pi { \alpha }^{ 2 } } \sqrt { \pi { \alpha }^{ 2 } } }{ \sqrt { \pi { \alpha }^{ 2 } } \pi { \alpha }^{ 2 } } exp\left( -\frac { { x_{ \tilde { m } } }^{ 2 } }{ { \alpha }^{ 2 } } \right) \\ =\frac { \pi { \alpha }^{ 2 } }{ \sqrt { \pi { \alpha }^{ 2 } } \pi { \alpha }^{ 2 } } exp\left( -\frac { { x_{ \tilde { m } } }^{ 2 } }{ { \alpha }^{ 2 } } \right) \\ =\frac { 1 }{ \sqrt { \pi { \alpha }^{ 2 } } } exp\left( -\frac { { x_{ \tilde { m } } }^{ 2 } }{ { \alpha }^{ 2 } } \right) \\ }$

となります。

GGX分布の場合は...すみません、自分の積分力の欠如により、Symbolab Math Solver - Step by Step calculator で自動で解くと、

www.symbolab.com

$\displaystyle { { P }^{ 2- }\left( x_{ \tilde { m } } \right) =\int _{ -\infty }^{ \infty }{ \frac { 1 }{ \pi { \alpha }^{ 2 }\left( 1+\frac { { x_{ \tilde { m } } }^{ 2 }+{ y_{ \tilde { m } } }^{ 2 } }{ { \alpha }^{ 2 } } \right) } dy_{ \tilde { m } } } \\ { P }^{ 2- }\left( x_{ \tilde { m } } \right) =\frac { \alpha ^{ 2 } }{ 2\left( x_{ \tilde { m } }^{ 2 }+\alpha ^{ 2 } \right) ^{ \frac { 3 }{ 2 } } } }$

のようになります。

Λ （Beckmann分布）

先のMicrofacet入門（１）では、もう一つおざなりにしていたのは $\Lambda$ です。

$\displaystyle { \Lambda \left( \mu \right) =\frac { 1 }{ \mu } \int _{ \mu }^{ \infty }{ P_{ 2 }\left( q \right) \left( q-\mu \right) dq } \\ \Lambda \left( \omega _{ o } \right) =\frac { 1 }{ \cot { \theta _{ o } } } \int _{ \cot { \theta _{ o } } }^{ \infty }{ P_{ 2 }\left( q \right) \left( q-\cot { \theta _{ o } } \right) dq } }$

こちらも、なんと、解析的に解くことができます。 Beckmann分布の場合、

$\displaystyle { { P }^{ 2- }\left( x_{ \tilde { m } } \right) =\frac { 1 }{ \sqrt { \pi { \alpha }^{ 2 } } } exp\left( -\frac { { x_{ \tilde { m } } }^{ 2 } }{ { \alpha }^{ 2 } } \right) \\ \Lambda =\frac { 1 }{ \cot { \theta _{ o } } } \int _{ \cot { \theta _{ o } } }^{ +\infty }{ \left( x_{ \tilde { m } }-\cot { \theta _{ o } } \right) { P }^{ 2- }\left( x_{ \tilde { m } } \right) dx_{ \tilde { m } } } \\ =\frac { 1 }{ \cot { \theta _{ o } } } \int _{ \cot { \theta _{ o } } }^{ +\infty }{ \left( x_{ \tilde { m } }-\cot { \theta _{ o } } \right) \frac { 1 }{ \sqrt { \pi { \alpha }^{ 2 } } } exp\left( -\frac { { x_{ \tilde { m } } }^{ 2 } }{ { \alpha }^{ 2 } } \right) dx_{ \tilde { m } } } \\ =\frac { 1 }{ \cot { \theta _{ o } } } \int _{ \cot { \theta _{ o } } }^{ +\infty }{ x_{ \tilde { m } }\frac { 1 }{ \sqrt { \pi { \alpha }^{ 2 } } } exp\left( -\frac { { x_{ \tilde { m } } }^{ 2 } }{ { \alpha }^{ 2 } } \right) -\cot { \theta _{ o } } \frac { 1 }{ \sqrt { \pi { \alpha }^{ 2 } } } exp\left( -\frac { { x_{ \tilde { m } } }^{ 2 } }{ { \alpha }^{ 2 } } \right) dx_{ \tilde { m } } } \\ =\frac { 1 }{ \cot { \theta _{ o } } } \int _{ \cot { \theta _{ o } } }^{ +\infty }{ x_{ \tilde { m } }\frac { 1 }{ \sqrt { \pi { \alpha }^{ 2 } } } exp\left( -\frac { { x_{ \tilde { m } } }^{ 2 } }{ { \alpha }^{ 2 } } \right) dx_{ \tilde { m } } } -\frac { 1 }{ \cot { \theta _{ o } } } \int _{ \cot { \theta _{ o } } }^{ +\infty }{ \cot { \theta _{ o } } \frac { 1 }{ \sqrt { \pi { \alpha }^{ 2 } } } exp\left( -\frac { { x_{ \tilde { m } } }^{ 2 } }{ { \alpha }^{ 2 } } \right) dx_{ \tilde { m } } } \\ =\frac { 1 }{ \cot { \theta _{ o } } } \int _{ \cot { \theta _{ o } } }^{ +\infty }{ x_{ \tilde { m } }\frac { 1 }{ \sqrt { \pi { \alpha }^{ 2 } } } exp\left( -\frac { { x_{ \tilde { m } } }^{ 2 } }{ { \alpha }^{ 2 } } \right) dx_{ \tilde { m } } } -\int _{ \cot { \theta _{ o } } }^{ +\infty }{ \frac { 1 }{ \sqrt { \pi { \alpha }^{ 2 } } } exp\left( -\frac { { x_{ \tilde { m } } }^{ 2 } }{ { \alpha }^{ 2 } } \right) dx_{ \tilde { m } } } \\ }$

この後、左と右の項をそれぞれまた Symbolab Math Solver - Step by Step calculatorの力を借りて、

$\displaystyle { a=\frac { 1 }{ \alpha \tan { \theta } } =\frac { \cot \theta }{ \alpha } }$

と置きながら...

左の項 www.symbolab.com

$\displaystyle { \frac { 1 }{ \cot { \theta _{ o } } } \int _{ \cot { \theta _{ o } } }^{ +\infty }{ x_{ \tilde { m } }\frac { 1 }{ \sqrt { \pi { \alpha }^{ 2 } } } exp\left( -\frac { { x_{ \tilde { m } } }^{ 2 } }{ { \alpha }^{ 2 } } \right) dx_{ \tilde { m } } } =\frac { \sqrt { \alpha ^{ 2 } } \tan \left( θ \right) e^{ -\frac { \cot ^{ 2 } \left( θ \right) }{ a^{ 2 } } } }{ 2\sqrt { \pi } } \\ \\ a=\frac { 1 }{ \alpha \tan { \theta } } =\frac { \cot \theta }{ \alpha } \\ \alpha =\frac { 1 }{ a\tan { \theta } } \\ \\ =\frac { \frac { 1 }{ a\tan { \theta } } \tan \left( θ \right) e^{ -{ a }^{ 2 } } }{ 2\sqrt { \pi } } \\ =\frac { 1 }{ 2a\sqrt { \pi } } e^{ -{ a }^{ 2 } } }$

右の項

www.symbolab.com

$\displaystyle { \int _{ \cot { \theta _{ o } } }^{ +\infty }{ \frac { 1 }{ \sqrt { \pi { \alpha }^{ 2 } } } exp\left( -\frac { { x_{ \tilde { m } } }^{ 2 } }{ { \alpha }^{ 2 } } \right) dx_{ \tilde { m } } } =\frac { 1 }{ 2 } -\frac { \alpha { erf }\left( \frac { \cot \left( θ \right) }{ \alpha } \right) }{ 2\sqrt { \alpha ^{ 2 } } } \\ \\ \\ a=\frac { 1 }{ \alpha \tan { \theta } } =\frac { \cot \theta }{ \alpha } \\ \\ =\frac { 1 }{ 2 } -\frac { \alpha { erf }\left( a \right) }{ 2\alpha } \\ =\frac { 1 }{ 2 } -\frac { { erf }\left( a \right) }{ 2 } \\ =\frac { 1-erf\left( a \right) }{ 2 } \\ }$

したがって、まとめると、

$\displaystyle { \Lambda \left( \omega _{ o } \right) =\frac { 1 }{ 2a\sqrt { \pi } } e^{ -{ a }^{ 2 } }-\frac { 1-erf\left( a \right) }{ 2 } \\ \Lambda \left( \omega _{ o } \right) =\frac { 1 }{ 2a\sqrt { \pi } } e^{ -{ a }^{ 2 } }+\frac { erf\left( a \right) -1 }{ 2 } }$

ただし、

$\displaystyle { a=\frac { 1 }{ \alpha \tan { \theta } } =\frac { \cot \theta }{ \alpha } }$

と求めることができました。

Λ （GGX分布）

GGX分布の場合のΛは、

ここでまたまたSymbolab Math Solver - Step by Step calculatorの力を借りて、

www.symbolab.com

$\displaystyle { { P }^{ 2- }\left( x_{ \tilde { m } } \right) =\frac { \alpha ^{ 2 } }{ 2\left( x_{ \tilde { m } }^{ 2 }+\alpha ^{ 2 } \right) ^{ \frac { 3 }{ 2 } } } \\ \Lambda =\frac { 1 }{ \cot { \theta _{ o } } } \int _{ \cot { \theta _{ o } } }^{ +\infty }{ \left( x_{ \tilde { m } }-\cot { \theta _{ o } } \right) \frac { \alpha ^{ 2 } }{ 2\left( x_{ \tilde { m } }^{ 2 }+\alpha ^{ 2 } \right) ^{ \frac { 3 }{ 2 } } } dx_{ \tilde { m } } } \\ \\ { P }^{ 2- }\left( x_{ \tilde { m } } \right) =\frac { \alpha ^{ 2 } }{ 2\left( x_{ \tilde { m } }^{ 2 }+\alpha ^{ 2 } \right) ^{ \frac { 3 }{ 2 } } } \\ \Lambda =\frac { 1 }{ \cot { \theta _{ o } } } \int _{ \cot { \theta _{ o } } }^{ +\infty }{ \left( x_{ \tilde { m } }-\cot { \theta _{ o } } \right) { P }^{ 2- }\left( x_{ \tilde { m } } \right) dx_{ \tilde { m } } } \\ =\frac { 1 }{ 2 } \tan \left( θ \right) \left( \sqrt { a^{ 2 }+\cot ^{ 2 } \left( θ \right) } -\cot \left( θ \right) \right) \\ =\frac { 1 }{ 2 } \tan \left( θ \right) \left( \sqrt { \left( \frac { 1 }{ \alpha \tan { \theta } } \right) ^{ 2 }+\left( \frac { 1 }{ \tan { \theta } } \right) ^{ 2 } } -\frac { 1 }{ \tan { \theta } } \right) \\ =\frac { 1 }{ 2 } \tan \left( θ \right) \left( \sqrt { \left( \frac { 1 }{ \tan { \theta } } \right) ^{ 2 }\left( \frac { 1 }{ \alpha } +1 \right) } -\frac { 1 }{ \tan { \theta } } \right) \\ =\frac { 1 }{ 2 } \tan \left( θ \right) \left( \frac { 1 }{ \tan { \theta } } \sqrt { \frac { 1 }{ \alpha } +1 } -\frac { 1 }{ \tan { \theta } } \right) \\ =\frac { 1 }{ 2 } \tan \left( θ \right) \left( \frac { 1 }{ \tan { \theta } } \left( \sqrt { \frac { 1 }{ \alpha } +1 } -1 \right) \right) \\ =\frac { 1 }{ 2 } \left( \sqrt { \frac { 1 }{ \alpha } +1 } -1 \right) \\ =\frac { -1+\sqrt { \frac { 1 }{ \alpha } +1 } }{ 2 } \\ }$

ただし、

$\displaystyle { a=\frac { 1 }{ \alpha \tan { \theta } } =\frac { \cot \theta }{ \alpha } }$

と求めることができました。

実装例

 inline double chi_plus(double x) {
        return x < 0.0 ? 0.0 : 1.0;
    }
    inline double D_Beckman(const Vec3 &n, const Vec3 &h, double alpha) {
        double cosTheta = glm::dot(n, h);
        if (cosTheta < 1.0e-9) {
            return 0.0;
        }
        double cosTheta2 = cosTheta * cosTheta;
        double cosTheta4 = cosTheta2 * cosTheta2;
        double alpha2 = alpha * alpha;
        double chi = chi_plus(cosTheta);

        // \tan { \theta  } =\pm \frac { \sqrt { 1-\cos { \theta  }  }  }{ \cos { \theta  }  } \\ \tan ^{ 2 }{ \theta  } =\frac { 1-\cos ^{ 2 }{ \theta  }  }{ \cos ^{ 2 }{ \theta  }  } 
        double tanTheta2 = (1.0 - cosTheta2) / cosTheta2;
        return chi * std::exp(-tanTheta2 / alpha2) / (glm::pi<double>() * alpha2 * cosTheta4);
    }
    inline double D_GGX(const Vec3 &n, const Vec3 &h, double alpha) {
        double cosTheta = glm::dot(n, h);
        if (cosTheta < 1.0e-9) {
            return 0.0;
        }
        double cosTheta2 = cosTheta * cosTheta;
        double cosTheta4 = cosTheta2 * cosTheta2;
        double alpha2 = alpha * alpha;
        double chi = chi_plus(cosTheta);

        // \tan { \theta  } =\pm \frac { \sqrt { 1-\cos { \theta  }  }  }{ \cos { \theta  }  } \\ \tan ^{ 2 }{ \theta  } =\frac { 1-\cos ^{ 2 }{ \theta  }  }{ \cos ^{ 2 }{ \theta  }  } 
        double tanTheta2 = (1.0 - cosTheta2) / cosTheta2;
        return chi / (glm::pi<double>() * alpha2 * cosTheta4 * Sqr(1.0 + tanTheta2 / alpha2));
    }

    inline double lambda_beckmann(double cosTheta, double alpha) {
        double tanThetaO = std::sqrt(1.0 - cosTheta * cosTheta) / cosTheta;
        double a = 1.0 / (alpha * tanThetaO);
        return (std::erf(a) - 1.0) * 0.5 + std::exp(-a * a) / (2.0 * a * std::sqrt(glm::pi<double>()));
    }
    inline double lambda_ggx(double cosTheta, double alpha) {
        double tanThetaO = std::sqrt(1.0 - cosTheta * cosTheta) / cosTheta;
        double a = 1.0 / (alpha * tanThetaO);
        return -1.0 + std::sqrt(1.0 + 1.0 / (a * a)) * 0.5;
    }

    inline double G_2_beckmann(const Vec3 &omega_i, const Vec3 &omega_o, const Vec3 &omega_h, const Vec3 &n, double alpha) {
        double numer = chi_plus(glm::dot(omega_o, omega_h)) * chi_plus(glm::dot(omega_i, omega_h));
        double denom = (1.0 + lambda_beckmann(glm::dot(omega_o, n), alpha) + lambda_beckmann(glm::dot(omega_i, n), alpha));
        return numer / denom;
    }

    inline double G_2_ggx(const Vec3 &omega_i, const Vec3 &omega_o, const Vec3 &omega_h, const Vec3 &n, double alpha) {
        double numer = chi_plus(glm::dot(omega_o, omega_h)) * chi_plus(glm::dot(omega_i, omega_h));
        double denom = (1.0 + lambda_ggx(glm::dot(omega_o, n), alpha) + lambda_ggx(glm::dot(omega_i, n), alpha));
        return numer / denom;
    }