透視投影変換行列を実装しようと思って色々調べると，いくつか式が出てきて，どれが正解なのか，どれかが間違っているのか，と悩んでしまうので，自分なりに納得できる形にまとめてみる．
透視投影変換とは
3D空間にあるオブジェクトを画面に表示する際に，人間の視覚と同じように，近くの物は大きく，遠くの物は小さく映るように座標を変換することである．
透視投影では，カメラに映る範囲を視錐台(View Frsutum)として表し，この範囲の座標$(x, y, z)$の各要素の値が-1から1の範囲に収まるような変換を行う．つまり，$-1 \leq x \leq 1$，$-1 \leq y \leq 1$，$-1 \leq z \leq 1$となる．なお，$z$については，$0 \leq z \leq 1$とする場合もある．
視錐台の指定方法
カメラから見える範囲を表す視錐台を指定する方法がいくつかあるが，ここでは視野(field of view)の水平角度$fov_x$と垂直角度$fov_y$，カメラから最前面までの距離$near$，最奥面までの距離$far$，画面の縦横比であるアスペクト比(縦を1とする)$aspect$を用いる．なお，$0 \leq fov_x \leq \pi$，$0 \leq fov_y \leq \pi$とする．
x座標の変換
まず，$x$座標の変換について考える．まず，水平視野角の半分を次のように定める．
\theta = \frac{fov_x}{2}\\
すると，ある$z$での$x$の最大値$x_{max}$は，$\tan\theta$を用いて次のように表せる．
\tan\theta = \frac{|x_{max}|}{|z|}\\
\therefore |x_{max}| = |z|\tan\theta
よって，ある$z$で視錐台の範囲内で$x$が取りうる範囲は，$-|z|\tan\theta \leq x \leq |z|\tan\theta$となるので，$|z|\tan\theta$で割ると，$-1 \leq \frac{x}{|z|\tan\theta} \leq 1$となる．つまり，$x$については$\frac{1}{|z|\tan\theta}$を掛ければ良いことが分かる．
y座標の変換
垂直視野角$fov_y$とアスペクト比(aspect ratio)は，画面の横幅$width$と高さ$height$を用いて，次のような関係になる．
\psi = \frac{fov_y}{2}\\
aspect = \frac{width}{height} = \frac{2|z|\tan\theta}{2|z|\tan\psi} = \frac{|z|\tan\theta}{|z|\tan\psi}\\
\therefore |z|\tan\psi = \frac{|z|\tan\theta}{aspect}
よって，視錐台の範囲内で$y$が取りうる範囲が$-|z|\tan\psi \leq y \leq |z|\tan\psi$より，$x$の場合と同じように，$\frac{1}{|z|\tan\psi}$を掛ければ良いことが分かる．これは，上記の式を使うと次のように$\theta$から求まるように変形できる．
\frac{1}{|z|\tan\psi} = \frac{1}{\frac{|z|\tan\theta}{aspect}} = \frac{1}{|z|\tan\theta}\cdot aspect
z座標の変換(右手座標系の場合)
$x$座標の変換でも，$y$座標の変換でも，$|z|$で割っていた．これは，座標を同次座標系で考え，$V = (x', y', z', w)$において，$w = |z|$とすれば実現できる．しかし，$z$座標を変換する場合，単純に係数を掛けるだけだと，自分自身で割ってしまうので，座標に関係なく常に係数と同じ値になってしまう．そこで，$z$は次のように変換することがポイントとなる．
z' = Az + B
この式を$|z|$で割り，$near$や$far$を代入した際の値から$A$や$B$を求める方程式を解く．
ここでは，右手座標系から右手座標系に変換する場合を考える．右手座標系の場合，カメラから見える範囲にあるオブジェクトの$z$座標は全て負になる．また，$z = -near$の場合，変換後の$z' = 1$となり，$z = -far$の場合に変換後の$z' = -1$となる．
\frac{Az + B}{|z|}
\frac{A\cdot(-near) + B}{near} = -A + \frac{B}{near} = 1
\frac{A\cdot(-far) + B}{far} = -A + \frac{B}{far} = -1
\left(-A + \frac{B}{near}\right) - \left(-A + \frac{B}{far}\right) = 1 - (-1)
\frac{B}{near} - \frac{B}{far} = 2
B\cdot\left(\frac{far - near}{far \cdot near}\right) = 2
B = \frac{2 \cdot far \cdot near}{far - near}
-A + \frac{\frac{2 \cdot far \cdot near}{far - near}}{near} = -A + \frac{2 \cdot far}{far - near} = 1
A = \frac{2 \cdot far}{far - near} - 1 = \frac{2 \cdot far - (far - near)}{far - near} = \frac{far + near}{far - near}
透視投影変換行列(右手座標系の場合)
ここまでの結果から，右手座標系において，ある同次座標系$V = (x, y, z, 1.0)$に対して右から掛ける透視投影変換行列は次のようになる．
M = \left(
\begin{array}{cccc}
\frac{1}{\tan\theta} & 0.0 & 0.0 & 0.0\\
0.0 & \frac{1}{\tan\theta}\cdot aspect & 0.0 & 0.0 \\
0.0 & 0.0 & \frac{far + near}{far - near} & -1.0 \\
0.0 & 0.0 & \frac{2 \cdot far \cdot near}{far - near} & 0.0
\end{array}
\right)
$w$が$z$の絶対値になるように，-1.0にすることがポイントになっている．
z座標の変換(左手座標系の場合)
左手座標系の場合，$x$座標や$y$座標に関する考え方は変わらない．そして，$z$座標に関しては，$z = near$で$z' = -1$となり，$z = far$で$z' = 1$となるので，$A$と$B$は次のように求まる．
\frac{A \cdot near + B}{near} = A + \frac{B}{near} = -1
\frac{A \cdot far + B}{far} = A + \frac{B}{far} = 1
\left(A + \frac{B}{far}\right) - \left(A + \frac{B}{near}\right) = 1 - (-1) = 2
\frac{B}{far} - \frac{B}{near} = 2
B \cdot \frac{near - far}{far \cdot near} = 2
B = \frac{2 \cdot far \cdot near}{near - far} = -\frac{2 \cdot far \cdot near}{far - near}
A + \frac{\frac{2 \cdot far \cdot near}{near - far}}{far} = A + \frac{2 \cdot near}{near - far} = 1
A = 1 - \frac{2 \cdot near}{near - far} = \frac{(near - far) - 2\cdot near}{near - far} = \frac{-near - far}{near - far} = \frac{far + near}{far - near}
透視投影変換行列(左手座標系の場合)
上記の結果から，左手座標系の場合の透視投影変換行列は次のようになる．
M = \left(
\begin{array}{cccc}
\frac{1}{\tan\theta} & 0.0 & 0.0 & 0.0\\
0.0 & \frac{1}{\tan\theta}\cdot aspect & 0.0 & 0.0 \\
0.0 & 0.0 & \frac{far + near}{far - near} & 1.0 \\
0.0 & 0.0 & -\frac{2 \cdot far \cdot near}{far - near} & 0.0
\end{array}
\right)
最後に
本当は図も載せた方が分かりやすいのだけれど，描くのが面倒．何か楽に描けるツールがあればよいのだけれど．