「カメラの較正」の版間の差分

2011年10月23日 (日) 20:36時点における最新版

平面を使った較正

まあ、定番の手法。撮影した物体の3次元的な形状を計算するだけなら何もやらなくていいけど、対象物の大きさも知るためには既知の大きさの物体(参照物体)を撮影してみてそれが画像上で何画素になるか調べる必要がある。

例えば一辺20 [cm]位の立方体を用意して各頂点が画像上のどこになるか観測すると、3次元空間でそれぞれの方向に20 [cm]動いたら画像上でどう変化するかがわかるので簡単に3次元座標と画像座標の関係を求められる。これが直感的で分かりやすいけど、その立方体を用意するのが大変。この立方体が長さと角度の基準になるから、精密に作らないと復元結果が歪んじゃう。もちろん画像上の点と3次元座標の関係が得られればいいので、立方体じゃ無くても、球とかタワシとかフィギュアとか3次元物体であれば何でも使えるけど、いずれにせよ精密に作るのは難しい。それに丁寧に扱わないといけないし。

そこで平面ですよ。奥さん。

真っ直ぐな板とプリンタさえあれば、誰でも簡単に作れちゃうスグレモノ。壊れたってすぐ作り直せるし、ぞんざいに扱っても大丈夫。投げようが踏もうが噛み付こうが問題なし。まあ、本当に精度良く較正しようとするなら、それなりに注意深く作る必要があるけどね。。

座標系の設定

座標系は自由に決められるから何でもいいけど、後々混乱の元になるのできちんと決めておくべし。ここでは、

右手座標系
回転は右ネジ
世界座標系＝カメラ座標系 (世界座標系の原点にカメラ)
画像座標系は右向きにx軸、下向きにy軸
カメラ座標系と画像座標系のy軸は平行・同じ向き
平面の法線(z軸)は平面の表側が正

とする (図参照)。

ホモグラフィ

i 回目の観測画像における j 番目の画像座標 $\tilde{\boldsymbol{m}}_{ij}$ と参照平面上の座標 $\tilde{\boldsymbol{p}}_{ij}$ はホモグラフィ $H_i$ を通して次の関係があるとする。

$\tilde{\boldsymbol{m}}_{ij} \propto H_i \tilde{\boldsymbol{p}}_{ij}$

$H_i$ は参照平面に設定した座標と観測した画像座標の組が4点以上あれば計算できる。あと、この時参照平面に設定した単位系が世界座標系の単位系になる。つまり、平面上の座標を適当な原点から[mm]で表現すれば世界座標系も[mm]になる。

カメラパラメータとの関係

世界座標系からカメラ座標系への変換行列を $T_C$ 、モデル座標系から世界座標系への変換行列を $T_M$ とすると、モデル座標系の3次元点 $\tilde{\boldsymbol{X}}$ は以下のように投影される。

$\tilde{\boldsymbol{m}} \propto (K | \boldsymbol{0}) T_C T_M \tilde{\boldsymbol{X}}$
ただし、
$K = \begin{bmatrix} \alpha & \gamma & u \\ 0 & \beta & v \\ 0 & 0 & 1 \end{bmatrix},\quad T_C = \begin{bmatrix} I & \boldsymbol{0} \\ 0 & 1 \end{bmatrix},\quad T_M = \begin{bmatrix} \boldsymbol{r}_1 & \boldsymbol{r}_2 & \boldsymbol{r}_3 & \boldsymbol{t} \\ 0 & 0 & 0 & 1 \end{bmatrix}$

で、平面の x 軸 (= $\boldsymbol{r}_{1i}$ )、 y 軸 (= $\boldsymbol{r}_{2i}$ )と世界座標系における平面の原点位置 (= $\boldsymbol{t}_i$ )から、平面上の点 $\tilde{\boldsymbol{p}}_{ij}$ は次のようになる。

$\begin{align} \tilde{\boldsymbol{m}_{ij}} &\propto K \begin{bmatrix} \boldsymbol{r}_{1i} & \boldsymbol{r}_{2i} & \boldsymbol{t}_i \end{bmatrix} \tilde{\boldsymbol{p}}_{ij} \\ &\propto H_i \tilde{\boldsymbol{p}}_{ij} \end{align}$

そんな訳で、カメラパラメータとホモグラフィは以下の関係が成り立つ。

$s_i \begin{bmatrix} \boldsymbol{r}_{1i} & \boldsymbol{r}_{2i} & \boldsymbol{t}_i \end{bmatrix} = K^{-1} H_i$
ただし、 $s_i$ は定数倍の係数。

左辺の $\boldsymbol{r}_{1i}, \boldsymbol{r}_{2i}$ は大きさ1で互いに直交するので、 $H_i = [\boldsymbol{h}_{1i} \ \boldsymbol{h}_{2i} \ \boldsymbol{h}_{3i}]$ とすると次を得る。

$\begin{cases} \boldsymbol{h}_{1i}^\text{T} \omega \boldsymbol{h}_{1i} = \boldsymbol{h}_{2i}^\text{T} \omega \boldsymbol{h}_{2i} = s_i^2 \\ \boldsymbol{h}_{1i}^\text{T} \omega \boldsymbol{h}_{2i} = 0 \end{cases}$
ここで、
$\omega = \left(K^{-1} \right)^\text{T} K^{-1}$
とした。

未知数は対称行列 $\omega$ の成分 (6自由度)のみで、1回の観測で得られる拘束式が2本なので、異なる位置・姿勢の参照平面を3回以上観測すれば解く事ができる。これによって $\omega$ が求まれば $K$ が計算できる。

$\omega$ の算出

$\begin{cases} \boldsymbol{h}_{1i}^\text{T} \omega \boldsymbol{h}_{1i} - \boldsymbol{h}_{2i}^\text{T} \omega \boldsymbol{h}_{2i} = 0 \\ \boldsymbol{h}_{1i}^\text{T} \omega \boldsymbol{h}_{2i} = 0 \end{cases}$

を $\omega$ について整理する。 $\boldsymbol{h}_{1i} = [h_{11i} \ h_{21i} \ h_{31i}]$ とすると、

$\begin{bmatrix} h_{11i}^2 - h_{12i}^2 & 2(h_{11i} h_{21i} - h_{12i} h_{22i}) & h_{21i}^2 - h_{22i}^2 & 2(h_{11i} h_{31i} - h_{12i} h_{32i}) & 2(h_{21i} h_{31i} - h_{22i} h_{32i}) & h_{31i}^2 - h_{32i}^2 \\ h_{11i} h_{12i} & h_{11i} h_{22i} + h_{21i} h_{12i} & h_{21i} h_{22i} & h_{11i} h_{32i} + h_{31i} h_{12i} & h_{21i} h_{32i} + h_{31i} h_{22i} & h_{31i} h_{32i} \end{bmatrix} \begin{bmatrix} \omega_{11} \\ \omega_{12} \\ \omega_{22} \\ \omega_{13} \\ \omega_{23} \\ \omega_{33} \end{bmatrix} = 0$

となるから、各観測に対して求めたホモグラフィ $H_i$ をこの式に従って並べて0固有値に対応する固有ベクトルを求めれば $\omega$ になっているはず。あ、この式は $\omega_{ij} = \omega_{ji}$ として並び替えただけなので、確かめてから使ってね。

外部パラメータの算出

$\omega$ から $K$ が計算できるので、

$\begin{bmatrix} \boldsymbol{r}_{1i} & \boldsymbol{r}_{2i} & \boldsymbol{t}_i \end{bmatrix} = \frac{1}{s_i} K^{-1} H_i$

から外部パラメータ( $\boldsymbol{r}_{1i}, \boldsymbol{r}_{2i}, \boldsymbol{r}_{3i} = \boldsymbol{r}_{1i} {\times} \boldsymbol{r}_{2i}, \boldsymbol{t}_i$ )が求まる。ここで $s_i$ は $|\boldsymbol{r}_{1i}| = |\boldsymbol{r}_{2i}| = 1$ を利用して定める。これは、各観測における平面の位置・姿勢を表す。

仕上げ

これで大体のパラメータが求まったけど、なんだかんだで誤差が大きかったりする。なので、求めたパラメータを初期値にして例えば再投影誤差を最小化したりする;

$e^2 = \sum_{ij} || \text{proj}(\boldsymbol{p}_{ij} | K, R_i, \boldsymbol{t}_i) - \boldsymbol{m}_{ij} ||^2 \rightarrow \min$
ただし $\text{proj}(\boldsymbol{p}_{ij} | K, R_i, \boldsymbol{t}_i)$ は平面上の点 $\boldsymbol{p}_{ij}$ を推定したパラメータで画像平面へ投影した座標とする。

歪みパラメータも推定するなら(適当な初期値を設定して)この式に入れとけばいいはず。

@@ 25行: / 25行: @@
 * カメラ座標系と画像座標系のy軸は平行・同じ向き
 * 平面の法線(z軸)は平面の表側が正
-とする (図参照のこと)。
+とする (図参照)。
 == ホモグラフィ ==
-i回目の観測画像におけるj番目の画像座標<math>\tilde{\boldsymbol{m}}_{ij}</math>と参照平面上の座標<math>\tilde{\boldsymbol{p}}_{ij}</math>はホモグラフィ<math>H_i</math>を通して次の関係があるとする。
+i 回目の観測画像における j 番目の画像座標<math>\tilde{\boldsymbol{m}}_{ij}</math>と参照平面上の座標<math>\tilde{\boldsymbol{p}}_{ij}</math>はホモグラフィ<math>H_i</math>を通して次の関係があるとする。
 <math>
@@ 65行: / 65行: @@
 </math>
-で、平面のx軸(= <math>\boldsymbol{r}_{1i}</math>)、y軸(= <math>\boldsymbol{r}_{2i}</math>)と世界座標系における平面の原点位置(= <math>\boldsymbol{t}_i</math>)から、平面上の点<math>\tilde{\boldsymbol{p}}_{ij}</math>は次のようになる。
+で、平面の x 軸 (= <math>\boldsymbol{r}_{1i}</math>)、 y 軸 (= <math>\boldsymbol{r}_{2i}</math>)と世界座標系における平面の原点位置 (= <math>\boldsymbol{t}_i</math>)から、平面上の点<math>\tilde{\boldsymbol{p}}_{ij}</math>は次のようになる。
 <math>
@@ 104行: / 104行: @@
 とした。
-未知数は対称行列<math>\omega</math>の成分(6自由度)のみで、1回の観測で得られる拘束式が2本なので、異なる位置・姿勢の参照平面を3回以上観測すれば解く事ができる。
+未知数は対称行列<math>\omega</math>の成分 (6自由度)のみで、1回の観測で得られる拘束式が2本なので、異なる位置・姿勢の参照平面を3回以上観測すれば解く事ができる。
 これによって<math>\omega</math>が求まれば<math>K</math>が計算できる。
@@ 155行: / 155行: @@
 <math>
-   e^2 = \sum_{ij} || \text{proj}(\boldsymbol{p}_{ij} | K, R_i, \boldsymbol{t}_i) - \boldsymbol{m}_{ij} ||^2
+   e^2 = \sum_{ij} || \text{proj}(\boldsymbol{p}_{ij} | K, R_i, \boldsymbol{t}_i) - \boldsymbol{m}_{ij} ||^2 \rightarrow \min
 </math><br />
 ただし<math>\text{proj}(\boldsymbol{p}_{ij} | K, R_i, \boldsymbol{t}_i)</math>は平面上の点<math>\boldsymbol{p}_{ij}</math>を推定したパラメータで画像平面へ投影した座標とする。
 歪みパラメータも推定するなら(適当な初期値を設定して)この式に入れとけばいいはず。