本文へスキップ

STUDIO SPOON

Web / App / SEO / Graphic

STUDIO SPOON

Web / App / SEO / Graphic

© STUDIO SPOON Inc.

Engineering

射影行列を求める

WebGPU3DCG
H.IFront-end Engineer

WebGPUやWebGL（Three.js）を使い3Dを表現したい場合に、必ず出てくる3DCGの概念としてカメラがあります。
カメラを定義することで、頂点を画面描画用の正規化デバイス座標（NDC）空間へ変換し、更に立体的に見せることができます。

実際にVertex Shaderでは、カメラはProjection MatrixとView Matrixという行列で表現され、頂点座標に作用させます。

pos_clip = projectionMatrix * viewMatrix * modelMatrix * vec4f(position, 1);

Projection Matrixは、射影行列と呼ばれます。
射影行列とは、カメラから見た3D空間を、画面に表示するための座標空間へ変換する行列です。

一般的な3DCGを扱う数学ライブラリには、射影行列を得る関数が用意されていると思いますが、「その行列がどのようにして求められたのか」については追い求めず、ずっと射影行列を使ってきました。

なんとなく使っていたこの魔法の行列を、今回は求めてみようと思います。

概要

射影行列には、正射影行列（Orthographic Projection Matrix）と透視投影行列（Perspective Projection Matrix）があります。

正射影と透視投影の違いは、奥行き方向へのスケールの違いの有無です。
正射影では、奥行きによるスケールの違いがなく、言い換えると実際の寸法や比率が保たれるため、工業図面やCADなどで使用されます。
透視投影では、奥行きによるスケールの違いがあるため、人間の目で見た景色と同じような現実に近い表現ができます。

左: 正射影, 右: 透視投影

この両方の射影行列について求めていきます。

正射影行列（Orthographic Projection Matrix）

正射影行列は、3D空間の指定した直方体領域を、画面描画用のNDC空間へ変換する行列です。

正射影の変換

x, y, zは以下のように表されます。

x∈[left,right]y∈[bottom,top]z∈[−near,−far]x \in [left, right] \\[0.5em] y \in [bottom, top] \\[0.5em] z \in [-near, -far] \\[0.5em]

これを、以下の正規化デバイス座標（NDC）の範囲へ変換する行列を作成します。

x∈[−1,1]y∈[−1,1]z∈[0,1]x \in [-1, 1] \\[0.5em] y \in [-1, 1] \\[0.5em] z \in [0, 1] \\[0.5em]
Z方向のNDC範囲

WebGPUでは、「0 ~ 1」の範囲をとります。
WebGLでは、「-1 ~ 1」の範囲をとります。

最終的な形状

正射影行列は、最終的に以下のように定義されます。

P=[2r−l00r+ll−r02t−b0t+bb−t001n−fnn−f0001]P = \begin{bmatrix} \frac{2}{r-l} & 0 & 0 & \frac{r+l}{l-r} \\[0.5em] 0 & \frac{2}{t-b} & 0 & \frac{t+b}{b-t} \\[0.5em] 0 & 0 & \frac{1}{n-f} & \frac{n}{n-f} \\[0.5em] 0 & 0 & 0 & 1 \\ \end{bmatrix} r:right,l:left,t:top,b:bottom,n:near,f:farr: right, \quad l: left, \quad t: top, \quad b: bottom, \quad n: near, \quad f: far

導出

正射影では、遠近による縮小がないため、NDCの成分はそれぞれ独立して以下のように考えることができます。

xndc=Axx+Bxyndc=Ayy+Byzndc=Azz+Bzx_{ndc} = A_{x}x + B_{x} \\[0.5em] y_{ndc} = A_{y}y + B_{y} \\[0.5em] z_{ndc} = A_{z}z + B_{z} \\[0.5em]

したがって、求める行列は、最初から以下の形であることがわかります。

P=[Ax00Bx0Ay0By00AzBz0001]P = \begin{bmatrix} A_{x} & 0 & 0 & B_{x} \\[0.5em] 0 & A_{y} & 0 & B_{y} \\[0.5em] 0 & 0 & A_{z} & B_{z} \\[0.5em] 0 & 0 & 0 & 1 \\ \end{bmatrix}

頂点を考えると、以下のようになります。

[Ax00Bx0Ay0By00AzBz0001][xyz1]\begin{bmatrix} A_{x} & 0 & 0 & B_{x} \\[0.5em] 0 & A_{y} & 0 & B_{y} \\[0.5em] 0 & 0 & A_{z} & B_{z} \\[0.5em] 0 & 0 & 0 & 1 \\ \end{bmatrix} \begin{bmatrix} x \\[0.5em] y \\[0.5em] z \\[0.5em] 1 \\ \end{bmatrix}

計算すると、

xclip=Axx+Bxyclip=Ayy+Byzclip=Azz+Bzwclip=1x_{clip} = A_{x}x + B_{x} \\[0.5em] y_{clip} = A_{y}y + B_{y} \\[0.5em] z_{clip} = A_{z}z + B_{z} \\[0.5em] w_{clip} = 1

wclip=1w_{clip} = 1なので、

xndc=xclipwclip=Axx+Bxyndc=yclipwclip=Ayy+Byzndc=zclipwclip=Azz+Bz\begin{align} x_{ndc} = \frac{x_{clip}}{w_{clip}} = A_{x}x + B_{x} \\[0.5em] y_{ndc} = \frac{y_{clip}}{w_{clip}} = A_{y}y + B_{y} \\[0.5em] z_{ndc} = \frac{z_{clip}}{w_{clip}} = A_{z}z + B_{z} \\[0.5em] \end{align}

1. X方向を求める

X方向では、xがleftの時に、

xndc=l→−1x_{ndc} = l → -1

xがrightの時に、

xndc=r→1x_{ndc} = r → 1

に変換したいです。この条件を(1)(1)式に入れると、以下の連立方程式が得られます。

Axl+Bx=−1Axr+Bx=1\begin{align} A_{x}l + B_{x} = -1 \tag{1.1} \\[0.5em] A_{x}r + B_{x} = 1 \tag{1.2} \\[0.5em] \end{align}

(1.2)−(1.1)(1.2) - (1.1)から、

Ax(r−l)=2Ax=2r−lA_{x}(r - l) = 2 \\[0.5em] A_{x} = \frac{2}{r - l}

(1.2)(1.2)にAxA_{x}を代入すると、

(2r−l)r+Bx=1Bx=1−2rr−l=r−l−2rr−l=r+ll−r\left(\frac{2}{r - l} \right)r + B_{x} = 1 \\[0.5em] B_{x} = 1 - \frac{2r}{r - l} = \frac{r - l - 2r}{r - l} = \frac{r + l}{l - r}

ここまでを行列で表すと、

P=[2r−l00r+ll−r0Ay0By00AzBz0001]P = \begin{bmatrix} \frac{2}{r - l} & 0 & 0 & \frac{r + l}{l - r} \\[0.5em] 0 & A_{y} & 0 & B_{y} \\[0.5em] 0 & 0 & A_{z} & B_{z} \\[0.5em] 0 & 0 & 0 & 1 \\ \end{bmatrix}

2. Y方向を求める

同様に、Y方向では、

xndc=b→−1xndc=t→1x_{ndc} = b → -1 \\[0.5em] x_{ndc} = t → 1 \\[0.5em]

に変換します。この条件を(2)(2)式に入れると、以下の連立方程式が得られます。

Ayb+By=−1Ayt+By=1\begin{align} A_{y}b + B_{y} = -1 \tag{2.1} \\[0.5em] A_{y}t + B_{y} = 1 \tag{2.2} \\[0.5em] \end{align}

これを解くと、

Ay=2t−bBy=t+bb−tA_{y} = \frac{2}{t - b} \\[0.5em] B_{y} = \frac{t + b}{b - t} \\[0.5em]

ここまでを行列で表すと、

P=[2r−l00r+ll−r02t−b0t+bb−t00AzBz0001]P = \begin{bmatrix} \frac{2}{r - l} & 0 & 0 & \frac{r + l}{l - r} \\[0.5em] 0 & \frac{2}{t - b} & 0 & \frac{t + b}{b - t} \\[0.5em] 0 & 0 & A_{z} & B_{z} \\[0.5em] 0 & 0 & 0 & 1 \\ \end{bmatrix}

3. Z方向を求める

同様に、Z方向では、

zndc=−n→0zndc=−f→1z_{ndc} = -n → 0 \\[0.5em] z_{ndc} = -f → 1 \\[0.5em]

に変換します。
ここで、WebGPUにおけるzndcz_{ndc}の範囲と、カメラの向きが通常-Zを向いていることに留意します。

zndc∈[0,1],z=−near,z=−farz_{ndc} \in [0, 1], \quad z = -near, \quad z = -far

この条件を(3)(3)式に入れると、以下の連立方程式が得られます。

−Azn+Bz=0−Azf+Bz=1\begin{align} -A_{z}n + B_{z} = 0 \tag{3.1} \\[0.5em] -A_{z}f + B_{z} = 1 \tag{3.2} \\[0.5em] \end{align}

これを解くと、

Az=1n−fBz=nn−fA_{z} = \frac{1}{n - f} \\[0.5em] B_{z} = \frac{n}{n - f} \\[0.5em]

行列で表すと、

[2r−l00r+ll−r02t−b0t+bb−t001n−fnn−f0001]\begin{bmatrix} \frac{2}{r-l} & 0 & 0 & \frac{r+l}{l-r} \\[0.5em] 0 & \frac{2}{t-b} & 0 & \frac{t+b}{b-t} \\[0.5em] 0 & 0 & \frac{1}{n-f} & \frac{n}{n-f} \\[0.5em] 0 & 0 & 0 & 1 \\ \end{bmatrix}

となり、最初に求めたかった行列PPと一致します。

コード上の記述

コード上では以下のように定義されます。

コード上の正射影行列
[
        2 / (r - l),                 0,           0, 0,
                  0,       2 / (t - b),           0, 0,
                  0,                 0, 1 / (n - f), 0,
  (r + l) / (l - r), (t + b) / (b - t), n / (n - f), 1,
]

WebGPU（WebGL）では、Column-major order（列優先）（新しいタブで開く）となるため、求めた行列PPが転置した形になります。

透視投影行列（Perspective Projection Matrix）

透視投影行列は、3D空間を人間の目やカメラに近い遠近感を持った形で、NDC空間へ変換する行列です。
正射影とは異なり、カメラから近い物体ほど大きく、遠い物体ほど小さく表示されます。

透視投影の変換

遠近感を作るためには、Z方向（つまり深度方向）の値で、X, Yを割ります。
また、Vertex Shaderの戻り値として渡すpositionはクリップ空間になっていて、Vertex Shader以降のパイプライン上でW成分で除算する（これを透視除算（Perspective Divide）という）ことによって、正規化デバイス座標空間（NDC）が得られます。
透視除算を利用する形で、遠近感を得るために割る値をW成分に設定すると都合が良くなります。

wclip=−zw_{clip} = -zなので、

xndc∝x−z=xwclipyndc∝y−z=ywclipx_{ndc} \propto \frac{x}{-z} = \frac{x}{w_{clip}} \\[.5em] y_{ndc} \propto \frac{y}{-z} = \frac{y}{w_{clip}} \\[.5em]

最終的な形状

透視投影行列は、最終的に以下のように定義されます。

P=[faspect0000f0000FN−FNFN−F00−10]P = \begin{bmatrix} \frac{f}{aspect} & 0 & 0 & 0 \\[.5em] 0 & f & 0 & 0 \\[.5em] 0 & 0 & \frac{F}{N - F} & \frac{NF}{N - F} \\[.5em] 0 & 0 & -1 & 0 \\ \end{bmatrix} f=tan⁡(π2−θ2),θ=fieldOfViewYInRadiansN:near,F:farf = \tan\left(\frac{\pi}{2} - \frac{θ}{2}\right), \quad θ = fieldOfViewYInRadians \\[.5em] N: near, \quad F: far

導出

以下は、透視投影図をYZ平面で見たときの図です。

透視投影の範囲図をYZ平面で見た時の図

1. Y方向を求める

上図から、「あるzzにおいて、画面の一番上にあたるyyはいくつか？」を考えます。

三角関数の関係から、

tan⁡(fovY2)=ymax−z\tan\left(\frac{fovY}{2}\right) = \frac{y_{max}}{-z} ymax=−z×tan⁡(fovY2)y_{max} = -z × \tan\left(\frac{fovY}{2}\right)

ffを以下のように定義すると、

f=1tan⁡(fovY2)f = \frac{1}{\tan\left(\frac{fovY}{2}\right)}

yyは、以下のように表せます。

ymax=−zfy_{max} = -\frac{z}{f}

yndcy_{ndc}は、[−1,1][-1, 1]の範囲になるので、任意のyyにおいてyndcy_{ndc}を考えると、

yndc=yymax=fy−zy_{ndc} = \frac{y}{y_{max}} = \frac{fy}{-z}

透視除算により、wclipw_{clip}は以下のように定義したので、

wclip=−zw_{clip} = -z

yclipy_{clip}は以下のように求まります。

yclip=yndc×wclip=fy−z×−z=fyy_{clip} = y_{ndc} \times w_{clip} = \frac{fy}{-z} \times -z \\[.5em] = fy

ここまでを行列で表すと、

P=[?0000f0000??00−10]P = \begin{bmatrix} ? & 0 & 0 & 0 \\[.5em] 0 & f & 0 & 0 \\[.5em] 0 & 0 & ? & ? \\[.5em] 0 & 0 & -1 & 0 \end{bmatrix}

クリップ座標を求める計算は以下のようになり、求めたyclipy_{clip}とwclipw_{clip}が行列PPで表されていることがわかります。

rclip=P[xyzw]=[?0000f0000??00−10][xyzw]r_{clip} = P \begin{bmatrix} x \\[.5em] y \\[.5em] z \\[.5em] w \end{bmatrix} = \begin{bmatrix} ? & 0 & 0 & 0 \\[.5em] 0 & f & 0 & 0 \\[.5em] 0 & 0 & ? & ? \\[.5em] 0 & 0 & -1 & 0 \end{bmatrix} \begin{bmatrix} x \\[.5em] y \\[.5em] z \\[.5em] w \end{bmatrix}

2. X方向を求める

画面のアスペクト比を以下のように置きます。

a=widthheighta = \frac{width}{height}

アスペクト比を用いると、任意のxxに対するxclipx_{clip}は、yclipy_{clip}と同様の手順で以下のように導出されます。

xmaxx_{max}についてアスペクト比を用いると、

ymax=−zfxmax=a×ymax=a−zfy_{max} = \frac{-z}{f} \\[.5em] x_{max} = a \times y_{max} = a \frac{-z}{f}

xndcx_{ndc}についてxmaxx_{max}を用いると、

xndc=xxmaxxndc=fx−azx_{ndc} = \frac{x}{x_{max}} \\[.5em] x_{ndc} = \frac{fx}{-az}

以上より、xclipx_{clip}は以下のように求まります。

wclip=−zxclip=xndc×wclip=fx−az×−zxclip=faxw_{clip} = -z \\[.5em] x_{clip} = x_{ndc} \times w_{clip} = \frac{fx}{-az} \times -z \\[.5em] x_{clip} = \frac{f}{a}x

ここまでを行列で表すと、

P=[fa0000f0000??00−10]P = \begin{bmatrix} \frac{f}{a} & 0 & 0 & 0 \\[.5em] 0 & f & 0 & 0 \\[.5em] 0 & 0 & ? & ? \\[.5em] 0 & 0 & -1 & 0 \end{bmatrix}

3. Z方向を求める

Z方向では、以下のような変換を考えます。

z=−near→zndc=0z=−far→zndc=1z = -near → z_{ndc} = 0 \\[.5em] z = -far → z_{ndc} = 1

行列では以下のように表されます。

[fa0000f0000AB00−10][xyz1]\begin{aligned} \begin{bmatrix} \frac{f}{a} & 0 & 0 & 0 \\[.5em] 0 & f & 0 & 0 \\[.5em] 0 & 0 & A & B \\[.5em] 0 & 0 & -1 & 0 \end{bmatrix} \begin{bmatrix} x \\[.5em] y \\[.5em] z \\[.5em] 1 \end{bmatrix} \end{aligned}

zclipz_{clip}は、以下のようになります。

zclip=Az+Bz_{clip} = Az + B

BBは任意の変数として、仮定されたものになります。

求めたいzndcz_{ndc}は、以下のような式で表されます。

zndc=zclipwclip=Az−z+B−z=−A−Bzs=−B,c=−Aと置くとzndc=sz+cz_{ndc} = \frac{z_{clip}}{w_{clip}} = \frac{Az}{-z} + \frac{B}{-z} = -A - \frac{B}{z} \\[1em] s = -B, \quad c = -A \quad \text{と置くと} \\[1em] z_{ndc} = \frac{s}{z} + c

この未知数ss, ccについて求めていきます。

nearnearでは、以下のように変換したかったので、

z=−near→zndc=0z = -near → z_{ndc} = 0

求める式は、

s−near+c=0(3.1)\frac{s}{-near} + c = 0 \tag{3.1}

farfarでは、以下のように変換したかったので、

z=−far→zndc=1z = -far → z_{ndc} = 1

求める式は、

s−far+c=1(3.2)\frac{s}{-far} + c = 1 \tag{3.2}

したがって、(3.1)(3.1)と(3.2)(3.2)の連立方程式を解けば、未知数ss, ccを求めることができます。

ここで、以下のように置きます。

near=N,far=Fnear = N, \quad far = F

ssについて、(3.2)−(3.1)(3.2) - (3.1)をすると、

s−F−s−N=1s×N−F×N−s×F−N×F=1sN−F−N×F=1s=NFF−N\frac{s}{-F} - \frac{s}{-N} = 1 \\[.5em] \frac{s \times N}{-F \times N} - \frac{s \times F}{-N \times F} = 1 \\[.5em] s\frac{N - F}{-N \times F} = 1 \\[.5em] s = \frac{NF}{F - N}

ccについて、(3.1)(3.1)にssを代入すると、

s−N+c=0c=s×1N=NFF−N×1Nc=FF−N\frac{s}{-N} + c = 0 \\[.5em] c = s \times \frac{1}{N} = \frac{NF}{F - N} \times \frac{1}{N} \\[.5em] c = \frac{F}{F - N}

以上より、zndcz_{ndc}は以下のようになります。

zndc=sz+czndc=NFF−N×1z+FF−Nz_{ndc} = \frac{s}{z} + c \\[.5em] z_{ndc} = \frac{NF}{F - N} \times \frac{1}{z} + \frac{F}{F - N}

zclipz_{clip}は、以下のように求まります。

zclip=zndc×wclip,wclip=−zzclip=(NFF−N×1z+FF−N)×−zzclip=FN−Fz+NFN−Fz_{clip} = z_{ndc} \times w_{clip}, \quad w_{clip} = -z \\[1em] z_{clip} = \left(\frac{NF}{F - N} \times \frac{1}{z} + \frac{F}{F - N}\right) \times -z \\[.5em] z_{clip} = \frac{F}{N - F}z + \frac{NF}{N - F}

行列に表すと、以下のようになり求めたかった透視投影行列と一致します。

P=[fa0000f0000FN−FNFN−F00−10]P = \begin{bmatrix} \frac{f}{a} & 0 & 0 & 0 \\[.5em] 0 & f & 0 & 0 \\[.5em] 0 & 0 & \frac{F}{N - F} & \frac{NF}{N - F} \\[.5em] 0 & 0 & -1 & 0 \\ \end{bmatrix}

コード上の記述

コード上では以下のように定義されます。

コード上の透視投影行列
const f = Math.tan(Math.PI * 0.5 - 0.5 * fieldOfViewYInRadians);
const rangeInv = 1 / (zNear - zFar);
 
[
  f / aspect, 0,                       0,  0,
           0, f,                       0,  0,
           0, 0,         zFar * rangeInv, -1,
           0, 0, zNear * zFar * rangeInv,  0,
]

カテゴリー

Engineering

タグ

この記事をシェアする
All Articles

セクションへのリンクRELATED ARTICLES

Engineering
WebGPU APIであると便利なライブラリ
  • WebGPU
Engineering
XDからFigmaへ、変換プラグイン「XD to Figma Converter」を作りました
  • Figmaプラグイン
  • デザイン
  • 実装
Engineering
Particle Simulationを最適化する
  • WebGPU
  • Compute Shader