OpenCVのカメラ内部パラメータをOpenGLのプロジェクション行列に変換する

OpenGLで描画した3Dの画像を用いてOpenCVで画像処理したかったのですが、カメラ内部パラメータとプロジェクション行列の相互変換の方法がわからなかったので、ここにメモしておきます。

内容についてですが、正直あまり自信がありません。
間違い等あれば指摘して下さると非常に助かります。

# 0. 結論
先に結論です。

OpenCVのキャリブレーションで求めたカメラ内部パラメータ行列 $K$ を

```math
K=\left[\begin{matrix}
f_{x} & 0 & c_{x}\\
0 & f_{y} & c_{y}\\
0 & 0 & 1
\end{matrix}\right]
```

とすると、OpenGLのプロジェクション行列 $P$ は

```math
P=
\left[\begin{matrix}
\frac{2 f_{x}}{w} & 0 & \frac{2 c_{x}}{w} - 1 & 0\\
0 & - \frac{2 f_{y}}{h} & - \frac{2 c_{y}}{h} + 1 & 0\\
0 & 0 & \frac{z_{f} + z_{n}}{z_{f} - z_{n}} & - \frac{2 z_{f} z_{n}}{z_{f} - z_{n}}\\
0 & 0 & 1 & 0
\end{matrix}\right]
```

となります。
$w, h$ は画像の解像度です。
$z_{n}, z_{f}$ は $z_{near}, z_{far}$ の略で、z軸で表示する範囲の最小値、最大値です。 $(0 < z_n < z_f)$

ここで注意ですが、このプロジェクション行列 $P$ はOpenCVの座標系に合わせたものです。
OpenCVとOpenGLではy軸とz軸の向き(符号)がそれぞれ逆です。
もし、OpenGLの座標系に合わせたプロジェクション行列を使用する場合、行列 $P$ の２列目と３列目の符号を全て反転します。

```math
P'=
P
\left[\begin{matrix}
1 & 0 & 0 & 0\\
0 & -1 & 0 & 0\\
0 & 0 & -1 & 0\\
0 & 0 & 0 & 1
\end{matrix}\right]
=
\left[\begin{matrix}
\frac{2 f_{x}}{w} & 0 & - \frac{2 c_{x}}{w} + 1 & 0\\
0 & \frac{2 f_{y}}{h} & \frac{2 c_{y}}{h} - 1 & 0\\
0 & 0 & - \frac{z_{f} + z_{n}}{z_{f} - z_{n}} & - \frac{2 z_{f} z_{n}}{z_{f} - z_{n}}\\
0 & 0 & -1 & 0
\end{matrix}\right]
```

また、プロジェクション行列 $P'$ を使用する場合は全ての頂点を次のように変換する必要があります。
OpenCVの点を $(X_{CV}, Y_{CV}, Z_{CV})$ 、OpenGLの点を $(X_{GL}, Y_{GL}, Z_{GL})$ とすると、

```math
\left[\begin{matrix}X_{GL}\\Y_{GL}\\Z_{GL}\end{matrix}\right]=
\left[\begin{matrix}
1 & 0 & 0\\
0 & -1 & 0\\
0 & 0 & -1\\
\end{matrix}\right]
\left[\begin{matrix}X_{CV}\\Y_{CV}\\Z_{CV}\end{matrix}\right]=
\left[\begin{matrix}X_{CV}\\-Y_{CV}\\-Z_{CV}\end{matrix}\right]
```

のように変換します。

以下はこの式の導出方法についての解説です。

# 1. OpenCVのカメラ内部パラメータについて
OpenCVを使ってカメラキャリブレーションを行うと、次のようなカメラ内部パラメータ行列 $K$ が求まります。

```math
K=\left[\begin{matrix}
f_{x} & 0 & c_{x}\\
0 & f_{y} & c_{y}\\
0 & 0 & 1
\end{matrix}\right]
```

ここで $f_x$ と $f_y$ はピクセル単位の焦点距離であり、レンズの中心からカメラセンサーまでの距離を示しています。（画素のアスペクト比が1:1ならば $f_x$ と $f_y$ は等しいので、通常は $f_x = f_y$ です。）

$c_x$ と $c_y$ はピクセル単位の光学中心座標です。通常は画像のおおよそ中央に位置します。

この行列 $K$ を用いることでカメラ座標系の任意の点 $(X, Y, Z)$ を画像座標系の点 $(u, v)$ に変換できます。

```math
s
\left[\begin{matrix}u\\v\\1\end{matrix}\right]=
\left[\begin{matrix}
f_{x} & 0 & c_{x}\\
0 & f_{y} & c_{y}\\
0 & 0 & 1
\end{matrix}\right]
\left[\begin{matrix}X\\Y\\Z\end{matrix}\right]
\tag{1}
```

ここで $s$　は右辺のz座標です。
今回の場合は $s=Z$ となります。
両辺をsで割り、行列を計算すると、

```math
\left[\begin{matrix}u\\v\\1\end{matrix}\right]=
\left[\begin{matrix}
f_{x} & 0 & c_{x}\\
0 & f_{y} & c_{y}\\
0 & 0 & 1
\end{matrix}\right]
\left[\begin{matrix}X/Z\\Y/Z\\1\end{matrix}\right]=
\left[\begin{matrix}
f_{x}\times\frac{X}{Z} + c_{x}\\
f_{y}\times\frac{Y}{Z} + c_{y}\\
1
\end{matrix}\right]
```

となります。

これが、OpenCVが `カメラ座標系の3次元の点` から `画像座標系の２次元の点` へと投影するモデルです。

# 2. OpenGLのプロジェクション行列へ変換
まずはOpenCVとOpenGLの描画範囲についてそれぞれ考えます。
描画する画面の解像度を $(w, h)$ とすると
OpenCVでは点 $(0, 0)$ が画面**左上**、点 $(w, h)$ が画面**右下**となるように描画されます。
OpenGLでは点 $(-1, -1)$ が画面**左下**、点 $(1, 1)$ が画面**右上**となるように描画されます。

そのため、OpenCVの画面上の点を $(u, v)$ 、OpenGLの画面上の点を $(u_{GL}, v_{GL})$ とすると、

```math
\left[\begin{matrix}u_{GL}\\v_{GL}\\1\end{matrix}\right]=
\left[\begin{matrix}
\frac{2}{w} & 0 & -1\\
0 & -\frac{2}{h} & 1\\
0 & 0 & 1
\end{matrix}\right]
\left[\begin{matrix}u\\v\\1\end{matrix}\right]
```

が成り立ちます。
これを式(1)と組み合わせると、

```math
s
\left[\begin{matrix}u_{GL}\\v_{GL}\\1\end{matrix}\right]=
\left[\begin{matrix}
\frac{2}{w} & 0 & -1\\
0 & -\frac{2}{h} & 1\\
0 & 0 & 1
\end{matrix}\right]
\left[\begin{matrix}
f_{x} & 0 & c_{x}\\
0 & f_{y} & c_{y}\\
0 & 0 & 1
\end{matrix}\right]
\left[\begin{matrix}X\\Y\\Z\end{matrix}\right]=
\left[\begin{matrix}
\frac{2 f_{x}}{w} & 0 & \frac{2 c_{x}}{w} - 1\\
0 & - \frac{2 f_{y}}{h} & - \frac{2 c_{y}}{h} + 1\\
0 & 0 & 1
\end{matrix}\right]
\left[\begin{matrix}X\\Y\\Z\end{matrix}\right]
```

となります。
また、OpenGLではz座標が -1 から 1 の範囲のみを描画します。( -1 が手前側です)
そのため、z座標を変換する必要があります。
$z_{n}, z_{f}$　( $0 < z_n < z_f$ ) をそれぞれz軸で表示する範囲の最小値、最大値とすると、

```math
t
\left[\begin{matrix}u_{GL}\\v_{GL}\\s\\1\end{matrix}\right]=
\left[\begin{matrix}
\frac{2 f_{x}}{w} & 0 & \frac{2 c_{x}}{w} - 1 & 0\\
0 & - \frac{2 f_{y}}{h} & - \frac{2 c_{y}}{h} + 1 & 0\\
0 & 0 & \frac{z_{f} + z_{n}}{z_{f} - z_{n}} & - \frac{2 z_{f} z_{n}}{z_{f} - z_{n}}\\
0 & 0 & 1 & 0
\end{matrix}\right]
\left[\begin{matrix}X\\Y\\Z\\1\end{matrix}\right]
```

となります。
ここで $t$　は右辺のw座標です。
今回の場合は $t=Z$ となります。
つまり、 $t$ は式(1)の $s$ と同じような役割を果たします。

以上より、OpenGLのプロジェクション行列 $P$ は

```math
P=
\left[\begin{matrix}
\frac{2 f_{x}}{w} & 0 & \frac{2 c_{x}}{w} - 1 & 0\\
0 & - \frac{2 f_{y}}{h} & - \frac{2 c_{y}}{h} + 1 & 0\\
0 & 0 & \frac{z_{f} + z_{n}}{z_{f} - z_{n}} & - \frac{2 z_{f} z_{n}}{z_{f} - z_{n}}\\
0 & 0 & 1 & 0
\end{matrix}\right]
```

となります。

**ただし、このプロジェクション行列 $P$ はOpenCVの座標系に合わせたものです。**
OpenGLの座標系に変換する場合については、その解説が以下に続きます。

# 3. OpenCVとOpenGLの座標系の違い
OpenCVとOpenGLは共にx軸のプラス方向が右向きです。
しかし、y軸とz軸の向き(符号)についてはOpenCVとOpenGLとでそれぞれ逆になります。
OpenCVではy軸のプラス方向が下向きであるのに対し、OpenGLではy軸のプラス方向が上向きです。
OpenCVではz軸のプラス方向が奥向きであるのに対し、OpenGLではz軸のプラス方向が手前向きです。
そのため、OpenGLの座標系に合わせたプロジェクション行列を使用する場合、行列 $P$ の２列目と３列目の符号を全て反転します。

```math
P'=
P
\left[\begin{matrix}
1 & 0 & 0 & 0\\
0 & -1 & 0 & 0\\
0 & 0 & -1 & 0\\
0 & 0 & 0 & 1
\end{matrix}\right]
=
\left[\begin{matrix}
\frac{2 f_{x}}{w} & 0 & - \frac{2 c_{x}}{w} + 1 & 0\\
0 & \frac{2 f_{y}}{h} & \frac{2 c_{y}}{h} - 1 & 0\\
0 & 0 & - \frac{z_{f} + z_{n}}{z_{f} - z_{n}} & - \frac{2 z_{f} z_{n}}{z_{f} - z_{n}}\\
0 & 0 & -1 & 0
\end{matrix}\right]
```

また、プロジェクション行列 $P'$ を使用する場合は全ての頂点を次のように変換する必要があります。
OpenCVの点を $(X_{CV}, Y_{CV}, Z_{CV})$ 、OpenGLの点を $(X_{GL}, Y_{GL}, Z_{GL})$ とすると、

```math
\left[\begin{matrix}X_{GL}\\Y_{GL}\\Z_{GL}\end{matrix}\right]=
\left[\begin{matrix}
1 & 0 & 0\\
0 & -1 & 0\\
0 & 0 & -1\\
\end{matrix}\right]
\left[\begin{matrix}X_{CV}\\Y_{CV}\\Z_{CV}\end{matrix}\right]=
\left[\begin{matrix}X_{CV}\\-Y_{CV}\\-Z_{CV}\end{matrix}\right]
```

のように変換します。

# 参考
カメラキャリブレーションと3次元再構成
http://opencv.jp/opencv-2.2/cpp/camera_calibration_and_3d_reconstruction.html

カメラ内部パラメータとは
https://mem-archive.com/2018/02/21/post-157/

OpenCVの内部パラメータでOpenGLの透視投影行列を作成（そのに）
https://13mzawa2.hateblo.jp/entry/2016/06/12/202640