Barycentric Coordinatesの計算とPerspective Correction, Partial Derivativeについて

Visibility Bufferなどのレンダリング技法を採用した際に必要になるかもしれない、Barycentric Coordinatesの計算と、それに関連する、Perspective CorrectionとPartial Derivativeについて考えてみたいと思います。

Visibility Bufferなどのレンダリング技法を採用した場合、従来のG-Bufferの代わりに、Render Targetには、Triangle IDなどの、ごくわずかな情報のみを格納します。そして、シェーディングを行う際は、格納された参照情報をもとに、実際にシェーディングを行うPrimitiveの様々な情報を読み出すと同時に、各種Attributeを補間するために、Primitive上での、シェーディングする点の重心座標:Barycentric Coordinatesを計算しなくてはなりません。また、Screen座標系における、各種AttributeのPartial Derivativeが必要になる場合もあります。GPUのGraphics Pipelineを使っているときは、何気なくRasterizerが計算してくれているのですが、そうでない場合は自分で計算する必要があります。

参照

Deferred Attribute Interpolation for Memory-Efficient Deferred Shading
http://cg.ivd.kit.edu/publications/2015/dais/DAIS.pdf

CLIPPING USING HOMOGENEOUS COORDINATES
https://www.microsoft.com/en-us/research/wp-content/uploads/1978/01/p245-blinn.pdf

Barycentric Coordinatesの計算

まず、レンダリングするPrimitiveの種類はいくつか考えられますし、重心座標系は三角形に限った話ではありませんが、ここでは三角形について考えます。三角形の各頂点をP1,P2,P3とすると、重心座標:(λ1,λ2,λ3)で表される位置:Pは、下記の様に計算できます。

\begin{aligned}    P = \lambda_1 P_1 + \lambda_2 P_2 + \lambda_3 P_3     \end{aligned}
(1)

また、三角形の平面上では、以下の式が成り立つので、重心座標の3変数のうちの2つが計算できれば十分です。

\begin{aligned}    \lambda_1 + \lambda_2 + \lambda_3 = 1    \end{aligned}
(2)

2次元の場合

まず初めに、2次元の三角形:(P1,P2,P3)と、その三角形上の点:Pが与えられた場合の、重心座標の計算方法について考えます。
P3を原点とした相対座標を以下のように定義します。

\begin{aligned}    D_1 =& P_1 - P_3 \\  D_2 =& P_2 - P_3 \\  D =& P - P_3 \\    \end{aligned}

そして、マトリクス:Tを以下の様に定義します。

\begin{aligned}    T =   \quad  \begin{pmatrix}   D_1.x & D_2.x \\  D_1.y & D_2.y   \end{pmatrix}  \quad    \end{aligned}

すると以下の式が導出できます。これは、初めに書いた重心座標の式(1),(2)の変形です。また、言い換えれば、D1,D2を基底とした、斜交座標系での座標を求めている式です。

\\  \begin{aligned}    T   \quad  \begin{pmatrix}   \lambda_1 \\  \lambda_2   \end{pmatrix}  \quad  =  \quad  \begin{pmatrix}   D.x \\  D.y   \end{pmatrix}  \quad    \end{aligned}

このことから、以下がの式が求まります。

\\  \begin{aligned}    \quad  \begin{pmatrix}   \lambda_1 \\  \lambda_2   \end{pmatrix}  \quad  =  T^{-1}  \quad  \begin{pmatrix}   D.x \\  D.y   \end{pmatrix}  \quad    \end{aligned}

つまり、Tの逆マトリクスが求まれば、重心座標系での位置が求まるというわけです。従って以下の様に計算します。

\\  \begin{aligned}    \det(T) =& D_1.x D_2.y - D_1.y D_2.x \\  \lambda_1 =& \frac{D.x D_2.y - D.y D_2.x}{\det(T)} \\  \lambda_2 =& \frac{D_1.x D.y - D_1.y D.x}{\det(T)} \\  \lambda_3 =& 1 - \lambda_1 - \lambda_2 \\    \end{aligned}

3次元の場合

2次元の場合と同様に、まず、P3を原点とした相対座標を以下のように定義します。

\\  \begin{aligned}    D_1 =& P_1 - P_3 \\  D_2 =& P_2 - P_3 \\  D =& P - P_3 \\    \end{aligned}

次に、これらを、D1とD2を基底とするアフィン平面に投影します。こうすることで、三角形が存在する平面上の2次元との座標として考えることができます。また、アフィン平面に投影しても、重心座標は保存されるようです。
具体的には以下の様に計算します。

\\  \begin{aligned}    D_{1x} =& dot(D_1, D_1) \\  D_{1y} =& D_{2x} = dot(D_1, D_2) \\  D_{2y} =& dot(D_2, D_2) \\  D_x =& dot(D, D_1) \\  D_y =& dot(D, D_2) \\    \end{aligned}

\\  \begin{aligned}    \det(T) =& D_{1x} D_{2y} - D_{1y} D_{2x} \\  \lambda_1 =& \frac{D_x D_{2y} - D_y D_{2x}}{\det(T)} \\  \lambda_2 =& \frac{D1_x D_y - D_1y D_x}{\det(T)} \\  \lambda_3 =& 1 - \lambda_1 - \lambda_2 \\    \end{aligned}

これで、3次元、2次元のいずれでも、重心座標を求めることができると思います。

Perspective Correctionについて

次に、View座標系の三角形をScreen座標系に透視変換した後の三角形における重心座標と、Attributeの補間について考えてみたいと思います。
まず、View座標系の三角形の各頂点と、重心座標を求めたい点をScreen座標に変換し、Screen座標系における重心座標を求めます。これは、既出の2次元における重心座標の求め方で計算できます。この時のスクリーン座標系における頂点と、重心座標をそれぞれ、Ps1,Ps2,Ps3, Ps:(λs1,λs2,λs3) とします。これらは明らかに、3次元のView座標系における頂点と重心座標、それぞれP1,P2,P3とP:(λ1,λ2,λ3)とは異なります。

まず、Screen座標に透視変換した三角形における重心位置の計算は以下の通りです。

\begin{aligned}    P_s = P_{s1}\lambda_{s1} + P_{s2}\lambda_{s2} + P_{s3}\lambda_{s3}    \end{aligned}

各頂点位置:Ps1,Ps2,Ps3そして、Psは、透視変換される際に、w成分(View座標系におけるz値)によって除算されているので、上式は以下と同じ意味を持ちます。(実際の計算では、View->Proj->Viewportと変換される過程でいろいろ計算が行われますが、非線形の演算は、w要素による除算のみです。)

\begin{aligned}    \frac{P}{P.w} = \frac{P_1}{P_1.w}\lambda_{s1} + \frac{P_2}{P_2.w}\lambda_{s2} + \frac{P_3}{P_3.w}\lambda_{s3}    \end{aligned}

つまり、Screen座標系における重心座標の計算は、3次元空間における各頂点位置の線形補間ではなく、各頂点位置をwで割ったものの線形補間と同義となります。この、View座標系における重心座標の計算とScreen座標系における重心座標の計算の関係は、頂点位置に限らず、すべてのAttribute、たとえばUVやColorなどにも当てはまります。
ここで、任意のAttribute:AをScreen座標系の重心位置で計算することを考えるならば、直接計算することは出来ませんが、A/P.wならば計算することは出来ます。

\begin{aligned}    \frac{A}{P.w} = \frac{A_1}{P_1.w}\lambda_{s1} + \frac{A_2}{P_2.w}\lambda_{s2} + \frac{A_3}{P_3.w}\lambda_{s3}    \end{aligned}

しかし、通常我々が実際に補間値として知りたい値は、A/P.wではなくAです。ここで、あるAttributeが、全ての頂点で1である場合を考えると、以下の式が求まります。

\begin{aligned}    \frac{1}{P.w} = \frac{1}{P_1.w}\lambda_{s1} + \frac{1}{P_2.w}\lambda_{s2} + \frac{1}{P_3.w}\lambda_{s3}    \end{aligned}

これを用いて、1/P.wを求めれば、A/P.wを求めた後に1/P.wで除算することで、求めたい値であるAが求まります。

Partial Derivativeについて

上記の、重心座標の計算と、Perspective Correctionによって、Screen座標系においても正しいAttributeの補間値が計算できると思います。しかし、これだけでは実際のシェーディングの際に不十分な場合があります。AttributeのPartial Derivativeは、たとえば、Textureをサンプリングするときに、u,v値のScreen座標系の各軸:x,yにおけるPartial Derivativeがわからないと、正しいフィルタリングを行うことができません。これらの値は、MipMapの選択や、Anisotropic Filteringのサンプリングを決定するときに必要です。
Partial Derivativeは、Rasterizerを使っているときは、隣接するPixelの計算結果との差分を計算することで、比較的容易に計算できますが、Visibility Bufferなどのレンダリング技法を採用しているときは、1Pixelごとに別々に計算する必要があります。

まずはじめに、Screen座標系における、重心座標についてもう一度考えてみます。ある三角形の透視変換後の各頂点位置とAttributeをPs1,Ps2,Ps3, A1,A2,A3とし、あるScreen座標:Psが、重心座標:(λs1,λs2,λs3)で表されるとき、Perspective Correctionを考慮しないAttributeの補間:Asは、以下の式で計算できます。

\begin{aligned}    P_s = P_{s1}\lambda_{s1} + P_{s2}\lambda_{s2} + P_{s3}\lambda_{s3} \\  A_s = A_1\lambda_{s1} + A_2\lambda_{s2} + A_3\lambda_{s3}    \end{aligned}

両式とも単純な線形補間の式なので、Asの各λにおけるPartial Derivativeは一定ですし、各λのスクリーン座標系のx,y各軸におけるPartial Derivativeの値も一定です。したがって、たとえば、Ps1を基準にした場合、上式のAsは以下のように書き換えられます。

\begin{aligned}    D_s =& P_s - P_{s1} \\  A_s =& A_1 \\  & + (D_s.x \frac{\partial \lambda_{s1}}{\partial x} + D_s.y\frac{\partial \lambda_{s1}}{\partial y}) A_1 \\  & + (D_s.x \frac{\partial \lambda_{s2}}{\partial x} + D_s.y\frac{\partial \lambda_{s2}}{\partial y}) A_2 \\  & + (D_s.x \frac{\partial \lambda_{s3}}{\partial x} + D_s.y\frac{\partial \lambda_{s3}}{\partial y}) A_3     \end{aligned}

さらに、以下のように書くことができます。

\begin{aligned}    A_s &= A_1 \\  & + D_s.x (\frac{\partial \lambda_{s1}}{\partial x}A_1 + \frac{\partial \lambda_{s2}}{\partial x}A_2 + \frac{\partial \lambda_{s3}}{\partial x}A_3) \\  & + D_s.y (\frac{\partial \lambda_{s1}}{\partial y}A_1 + \frac{\partial \lambda_{s2}}{\partial y}A_2 + \frac{\partial \lambda_{s3}}{\partial y}A_3)     \end{aligned}

ここからは簡便のためx軸についてのみ書きますがy軸についても同様です。

\begin{aligned}    \frac{\partial A_s}{\partial x} = \sum_{i=1}^{3}{\frac{\partial \lambda_{si}}{\partial x}A_i}    \end{aligned}

Perspective Correctionを考慮するならば、以下の式が考えられます。

\begin{aligned}    \frac{\partial \frac{A}{P.w}}{\partial x} =& \sum_{i=1}^{3}{\frac{\partial \lambda_{si}}{\partial x}\frac{A_i}{P_i.w}} \\  \frac{\partial \frac{1}{P.w}}{\partial x} =& \sum_{i=1}^{3}{\frac{\partial \lambda_{si}}{\partial x}\frac{1}{P_i.w}}    \end{aligned}

つまり、重心座標:(λs1,λs2,λs3)のScreen座標系の各軸:x,yにおけるPartial Derivativeの値が求まれば、AttributeのPartial Derivativeが計算できると思います。これらは、2次元における重心座標の計算の式から容易に求まります。少し繰り返しになりますが、以下の通りです。

\begin{aligned}    D_{s1} =& P_{s1} - P_{s3} \\  D_{s2} =& P_{s2} - P_{s3} \\  D_s =& P_s - P_{s3}    \end{aligned}

\begin{aligned}    \det(T) =& D_{s1}.x D_{s2}.y - D_{s1}.y D_{s2}.x \\    \frac{\partial \lambda_{s1}}{\partial x} =& \frac{D_{s2}.y}{\det(T)} \\  \frac{\partial \lambda_{s1}}{\partial y} =& \frac{-D_{s2}.x}{\det(T)} \\    \frac{\partial \lambda_{s2}}{\partial x} =& \frac{-D_{s1}.y}{\det(T)} \\  \frac{\partial \lambda_{s2}}{\partial y} =& \frac{D_{s1}.x}{\det(T)} \\    \frac{\partial \lambda_{s3}}{\partial x} =& \frac{P_{s1}.y-P_{s2}.y}{\det(T)} \\  \frac{\partial \lambda_{s3}}{\partial y} =& \frac{P_{s2}.x-P_{s1}.x}{\det(T)}    \end{aligned}

これで、Screen座標系において、Perspective Correctionを考慮したAttributeの補間を、Partial Derivativeを含めて計算できるようになりました。しかし、これは、三角形がView座標系からScreen座標系への変換の過程で、Near Clippingされることを考慮していません。つまり、三角形のいずれかの頂点のw成分が負になってしまう場合はうまく計算ができません。
これについては、参照論文の4.2章で軽く触れられていますが、詳細は記述されていません。(we use a simplified algorithm and only shrink the triangle to fit inside the clipping volume instead of creating new triangles.)

Near Clipが発生した場合について少し考えてみます。Screen座標系でシェーディングしなくてはならない点があるということは、三角形の少なくとも1つの頂点は、wが正であると言えると思います。この点を基準(Ps3)にしてλs1,λs2のPartial Derivativeを計算します。Ps1,Ps2に関してはNear Clippingが適用されているかもしれないと考えます。Near Clippingが適用された場合は、Ps1,Ps2についてView座標系でClippingされた新しい位置を計算する必要があります。加えて、その位置でのAttributeの値を計算します。これはView座標系なので、単純な線形補間で計算できます。
そして、更新した頂点を用いて、λs1,λs2を求めるために先の逆マトリクスを計算する方法を使った場合、λs1,λs2の合計が1を超える場合があります。これは当然で、ひとつの頂点のみNear Clippingした三角形は、(Near Clippingnのみを考慮した場合)Screen座標系においては四角形以上の頂点数を持つためです。
しかし、依然として計算されたλs1,λs2は、2本の基底ベクトル:Ds1,Ds2によって示される、Ps3からの相対座標であることに変わりはありません。加えて、λs3は負の値になる場合が考えられますが、三角形が存在する平面上に点がある限り、(2)の式で計算できます。(この場合は重心座標とは呼ぶことができないかもしれませんが、平面の方程式ではあります。)
したがって、各λsを求める式に変更がないので、Partial Derivativeを求める式にも変更がありません。
しかし、Near Clippingを考慮して、Partial Derivativeを求めることは、非常にコストが高いことがわかります。

簡易なPartial Derivativeの求め方(1)

Partial Derivativeをx,y各軸について正しく求めることをあきらめて、TextureのMipMapの選択だけでも簡単に解決できないでしょうか。
View座標系で、投影平面に正対するz=1の位置にある平面上に投影されたPixelサイズ相当の円の直径をdとした場合、任意の距離の正対する平面上でのPixelサイズ相当の円の直径は、d*P.zで計算できると思います。投影される三角形の平面が投影平面に対して傾けば、楕円状に投影されます。このときの三角形の法線をNとして、ViewベクトルをVとすれば、距離の変化による歪みを考えなければ(つまり局地的に正射影と考えた場合)、以下の式で楕円の長径が計算できると思います。

\begin{aligned}    d_{proj} = \frac{d P.z}{dot(N,V)}    \end{aligned}

また、三角形平面上での、楕円の長径のベクトルは、Viewベクトルが三角形平面上に投影されたものと同じ方向のはずなので、楕円の長径のベクトルは以下の様に計算できると思います。

\\  \begin{aligned}    \vec{d_{proj}} = \frac{d P.z}{dot(N, V)} normalize(V - N dot(N,V))    \end{aligned}

次に、このベクトルを、三角形の重心座標系に投影します。計算は既出の式と同様です。三角形の頂点から形成されるベクトル、D1とD2を基底とするアフィン平面に投影します。

\\  \begin{aligned}    D_1 =& P_1 - P_3 \\  D_2 =& P_2 - P_3 \\  D =& \frac{R P.w}{dot(N, v)} normalize(V - N dot(N,V))    \end{aligned}

\\  \begin{aligned}    D_{1x} =& dot(D_1, D_1) \\  D_{1y} =& D_{2x} = dot(D_1, D_2) \\  D_{2y} =& dot(D_2, D_2) \\  D_x =& dot(D, D_1) \\  D_y =& dot(D, D_2) \\    \end{aligned}

\\  \begin{aligned}    \det(T) =& D_{1x} D_{2y} - D_{1y} D_{2x} \\  \lambda_1 =& \frac{D_x D_{2y} - D_y D_{2x}}{\det(T)} \\  \lambda_2 =& \frac{D1_x D_y - D_1y D_x}{\det(T)} \\  \lambda_3 =& 1 - \lambda_1 - \lambda_2 \\    \end{aligned}

これは、View座標系における重心座標なので、そのまま任意のAttributeと乗算すれば、楕円の長径に相当するベクトルで、Attributeがどのように変化するかを計算することができると思います。これをPixel空間におけるPartial Derivativeとして代用すれば、Anisotropic Fliteringは出来ませんが、TextureのMipMapぐらいは決められるかと思います。また、これらの計算はView座標系で行うので、Near Clippingの心配をする必要がありません。

簡易なPartial Derivativeの求め方(2)

もう少しコストをかけてもいいケースを考えてみたいと思います。
Screen座標系におおけるPixel:(xp, yp, 0)について計算しているとします。このPixelから x方向に1Pixel分となりの座標は:(xp+1, yp, 0)で表すことができると思います。これらをView座標系に逆変換すれば、それぞれの位置が、View座標系のNear Clipping平面上で求まると思います。
加えて、これらの差分のベクトルを計算できると思います。この差分ベクトルを、(xv, yv, 0)とします。このベクトルの、透視変換後の見かけの大きさを変えず、シェーディングを行っている点の深度まで移動すれば、Dx:(xv*P.z/nz, yv*P.z/nz, 0)となるはずです。
このベクトルをView座標系の三角形平面上に、Viewベクトルに沿って投影します。これは、三角形の平面と、始点が(P+Dx)で向きがViewベクトル:Vの直線との交差位置を求めるのと同じです。

\\  \begin{aligned}  D_{px} = P+D_x + V\frac{-dot(D_x, N)}{dot(N, V)}    \end{aligned}

あとは、このDpxのPartial Derivativeを簡易なPartial Derivativeの求め方(1)で行った計算と同様に行えば良いはずです。同様にy軸方向も計算できるはずです。これで、Screen座標系におけるx,y軸のPartial Derivativeをそれぞれ計算することができます。
Screen座標系でPerspective Correction考慮に入れた計算よりは処理は少ないと思います。またNear Clippingを考慮にいれなくても問題ないと思われます。
本来ならば、隣のPixelのViewベクトルは異なるので、この計算は正確とは言えません。しかし、ベクトルの向きを変えると、三角形の平面と交差しなかったり、始点から負の方向で交差する場合があります。そのため、この計算を用いた方が安定的に値を求めることができると思います。

広告