2 3 2 3 2 3 Latest Algorithm for 3-D Reconstruction from Two Views Kento Yamada, Yasushi Kanazawa, Kenichi Kanatani 2 and Yasuyuki Sugaya 3 This paper presents a new algorithm for reconstructing the 3-D shape of the scene from point correspondences over two views. The basic principle is well known, but we incorporate into it the latest results of the authors studies of statistical optimization techniques. As a result, not only the accuracy increases but also the organization of the program becomes simpler, allowing this system to be used in a wider range of practical application. Also, our system take into consideration the degenerate situation where the camera is in a fixating configuration. This makes our system very practical, since this is the situation most frequently encountered in practice. Finally, we demonstrate the performance of our system using real images. Department of Knowledge-based Information Engineering, Toyohashi University of Technolgy 2 Department of omputer Science, Okayama University 3 Department of Information and omputer Sciences, Toyohashi University of Technology. 2 3 98, 9 (structure from motion) 5),7) 3 5) 2 3 9) ( ) 9) FNS ) ) ( 2 ) 9) Taubin 8) 8) ( 3 ) 9) 2 Bougnoux ) 9) ( 4 ) 9) 3 3) 3 2. : (x α, y α ), (x α, y α), α =,..., N ( 8). c 29 Information Processing Society of Japan
: 3 (X α, Y α, Z α), α =,..., N (, ) xy x y 2 5) 2 3 x y z x y xyz 2 3. (x α, y α), (x α, y α) F = (F ij) 9 u u = (F, F 2, F 3, F 2, F 22, F 23, F 3, F 32, F 33 ) () 2 3. 2 9) Taubin 8) 8) N [ ] N [a] = a/ a 8),9) 2 x y 3 Taubin 2 FNS ( ) (x α, y α), (x α, y α) 8 z α 3 z α = (x αx α, x αy α, f x α, y αx α, y αy α, f y α, f x α, f y α) (2) ( 2 ) 8 8 V [z α] x 2 α + x 2 α x αy α f x α x αy α f x α x αy α x 2 α + y α 2 f y α x αy α f x α f x α f y α f 2 x αy α yα 2 + x 2 α x αy α f x α f y α V [z α ] = x αy α x αy α yα 2 + y α 2 f y α f y α B f x α f y α f 2 @ f x α f y α f 2 A f x α f y α f 2 ( 3 ) 8 z, z α z = N z α, z α = z α z (4) N ( 4 ) 8 8 M TB, L TB N N M TB = z α z α, L TB = V [z α ] (5) ( 5 ) λ v M TB v = λl TB v (6) 3 f 2) f = 6 (3) 2 c 29 Information Processing Society of Japan
( 6 ) 9 ( u) u = N [ v (v, z)/f 2 ] (7) 3.2 u (x α, y α), (x α, y α) (ˆx α, ŷ α), (ˆx α, ŷ α) 3) 2) ( ) E = ˆx α = x α, ŷ α = y α, ˆx α = x α, ŷ α = y α, x α = ỹ α = x α = ỹ α =. (8) ( 2 ) 9 ξ α 9 9 V [ξ α ] ˆx α ˆx α + ˆx α x α + ˆx α x α ˆx αŷ α + ŷ α x α + ˆx αỹ α f (ˆx α + x α) ŷ αˆx α + ˆx αỹ α + ŷ α x α ξ α = ŷ αŷ α + ŷ αỹ α + ŷ αỹ α f (ŷ α + ỹ α) B f (ˆx α + x α) @ f (ŷ α + ỹ α) A V [ξ α ] = f 2 ˆx 2 α + ˆx 2 α ˆx αŷ α f ˆx α ˆx α ŷ α f ˆx α ˆx αŷ α ˆx 2 α + ŷ α 2 f ŷ α ˆx α ŷ α f ˆx α f ˆx α f ŷ α f 2 ˆx αŷ α ŷα 2 + ˆx 2 α ˆx αŷ α f ˆx α f ŷ α ˆx αŷ α ˆx αŷ α ŷα 2 + ŷ α 2 f ŷ α f ŷ α f ˆx α f ŷ α f 2 B f ˆx α f ŷ α f 2 @ f ˆx α f ŷ α f 2 A ( 3 ) EFNS 3.3 u ( 4 ) x α, ỹ α, x α, ỹ α! x α (u, ξ α) u u 2 u 3 ỹ α (u, V [ξ α ]u) u 4 u 5 u 6! B @ ˆx α ŷ α f A, (9) () x α ỹ α! (u, ξ α) (u, V [ξ α ]u) u u 4 u 7 u 2 u 5 u 8 ( 5 ) ˆx α, ŷ α, ˆx α, ŷ α! ˆx α B @ ŷ α A () f ˆx α x α x α, ŷ α y α ỹ α, ˆx α x α x α, ŷ α y α ỹ α (2) ( 6 ) E E = N ( x 2 α + ỹα 2 + x α N 7 2 + ỹ α 2 ) (3) ( 7 ) E E E E < 4 F u u 2 u 3 B F = @ u 4 u 5 u 6 A (4) u 7 u 8 u 9 E E 2 3.3 4 EFNS EFNS ) ( ) 9 9 M, L N ξ M = α ξ N α (u, V [ξ α ]u), L = (u, ξ α ) 2 V [ξ α ] (u, V [ξ α ]u) 2. (5) ( 2 ) 9 u 9 9 P u u 5 u 9 u 8 u 6 u 6 u 7 u 9 u 4 u 4 u 8 u 7 u 5 u 8 u 3 u 2 u 9 u = N [ u 9 u u 3 u 7 ], u 7 u 2 u u 8 B u 2 u 6 u 5 u 3 @ u 3 u 4 u 6 u A P u = I u u. (6) u u 5 u 4 u 2 ( 3 ) 9 9 X, Y X = M L, Y = P u XP u. (7) 7) 3 c 29 Information Processing Society of Japan
( 4 ) Y 9 v, v 2 ( 5 ) 9 û û = (u, v )v + (u, v 2)v 2 (8) ( 6 ) 9 u u = N [P u û]. (9) ( 7 ) u u u u < 6 u u u N [u + u ] EFNS 2) 4. F 2 f, f 9) f, f f = f 4. F F, F F 3 e, e 2 f, f k = (,, ) ξ = F k 2 (k, F F F k) e k 2 /(k, F k) e k 2 F k 2 (k, F k) 2, η = F k 2 (k, F F F k) e k 2 /(k, F k) e k 2 F k 2 (k, F k) 2 (2) f = f, f = f (2) + ξ + η ) 2 e, e (f e /e 3, f e 2/e 3) 2 (f e /e 3, fe 2 /e 3 ) 2 X f y f y O O 3 Y x Z 2 (2) (2) Bougnoux ) 8) (2) (F 33 =) (k, F k) = k (k, F k) = 2 2 3 8),9) (k, F k) 3 4.2 2 ξ η ξ = η = (H + H 2 )ξ + (H 22 + H 2 )η. (22) H + 2H 2 + H 22 H ij ( ) 2, H = 2(k, F k) 4 η 2 + 4(k, F k) 2 F k 2 η + 2 F k 4 (k, F k) 2 η + F k 2 ( ) 2, H 22 = 2(k, F k) 4 ξ 2 + 4(k, F k) 2 F k 2 ξ + 2 F k 4 (k, F k) 2 ξ + F k 2 ) H 2 = 4(k, F k) 4 ξη + 4(k, F k) ( F 2 k 2 ξ + F k 2 η + 4(k, F k)(k, F F F k) 3 (k, F k) <. min( F k, F k )/f Z x X Y 4 c 29 Information Processing Society of Japan
) ((k, F k) 2 ξ + F k )((k, 2 F k) 2 η + F k 2 (k, F k) 2 ( (k, F k) 2 ξη + F k 2 ξ + F k 2 η + F 2 ) (23) Z Z (22) ξ, η (2) f, f 4.3 (k, F k) = f = f 8),9) ( ) a a 5 a = (k, F k)4 2 a 2 = (k, F k) 2 ( F k 2 + F k 2 ) a 3 = 2 ( F k 2 F k 2 ) 2 + (k, F k)(4(k, F F F k) (k, F k) F 2 ) a 4 = 2( F F k 2 + F F k 2 ) ( F k 2 + F k 2 ) F 2 a 5 = F F 2 2 F 4 (24) ( 2 ) K(ξ) K(ξ) = a ξ 4 + a 2ξ 3 + a 3ξ 2 + a 4ξ + a 5 (25) ( 3 ) (k, F k) ξ ξ = a4 2a 3 (26) ( 4 ) 3 K (ξ) = 4a ξ 3 + 3a 2 ξ 2 + 2a 3 ξ + a 4 = (27) ( 5 ) (27) ξ 3 ξ 3 ξ 2 ξ ξ ξ 3 K(ξ 3 ) < K(ξ ) K(ξ 3 ) ξ = ξ 3 K(ξ 3) < K(ξ ) ( 6 ) f, f f = f = (28) f + ξ (29) t t (a) (b) 4 (a) (b) 4.4 f = f 2 4(a) 4(b) 8) 3 (2), (29) f, f 9) 3 3) 5 c 29 Information Processing Society of Japan
5. F f, f 2 t R {t, R} t t = 98 6),7) ( ) E E = diag(,, f f )F diag(,, f f ) (3) ( 2 ) EE t ( 3 ) (x α, y α ), (x α, y α) x α/f x α/f B x α = @ y α/f A, x B α = @ y α/f A (3) ( 4 ) t N t, x α, Ex α < (32) a, b, c a, b, c ( 5 ) t E t E = U diag(σ, σ 2, σ 3 )V. (33) t E t E t E ( 6 ) R R = U diag(,, det(uv ))V (34) t 2 4 E 6),7) 6. 3 6. {t, R} 3 5 (a) (b) (a) (b) 5(a) 5(b) 4),7) 3) {t, R} F F = diag(,, f f )(t R)diag(,, f f ) (35) () 9 3.2 3 u 7 F (ˆx α, ŷ α ), (ˆx α, ŷ α), α =,..., N {t, R} f, f (ˆx α, ŷ α ), (ˆx α, ŷ α) E 3.2 f = f f = f 3.2 E E F f = f f = f E 6.2 3 (ˆx α, ŷ α ), (ˆx α, ŷ α) 3 (X α, Y α, Z α ) ( ) 2 3 4 P, P P = diag(,, f ( ) f ) I, P = diag(,, f ( f ) R ) R t f f E (36) 6 c 29 Information Processing Society of Japan
情報処理学会研究報告 次の連立 次方程式を解いて 3 次元位置 (Xα, Yα, Zα ), α =,..., N を計算する (2) Pij, Pij はそれぞれ P, P の (ij) 要素 x α P3 f P B B y α P3 f P2 B @ x α P3 f P P f P y α 2 3 PN (3) x α P32 f P2 y α P32 f P22 f P x α P32 2 P f P y α 22 32 x α P33 f P3 y α P33 f P23 f P x α P33 3 P f P y α 23 33 x α P34 f P4 Xα B B B y α P34 f P24 @ Yα A = B A @ x α P34 f P4 Zα P f P y α 24 34 (37) A sgn(zα ) < なら sgn(x) は符号関数であり x >, x =, x < に応じ 図6 実画像シーンから抽出した対応点とそのフロー てそれぞれ,, を返す t とすべての (Xα, Yα, Zα ) の符号を換える 式 (36) は第 カメラのレンズ中心を原点とし その光軸を Z 軸とする XY Z 座標系を 仮定するものであり 復元された点の座標 (Xα, Yα, Zα ) はこの座標系に関する記述となる ステップ 3 は鏡像解 カメラの後方にある反転した形状 を除去するものである 式 (32) の判定で t の符号を E に適合させているが 基本行列 E は式 (3) により基礎行列 F から 定まる しかし 基礎行列 F は符号が不定であり この符号が正しくないと鏡像解が生じ る そこで全特徴点が第 カメラの前方にあるように t と 3 次元位置 (Xα, Yα, Zα ) の符号 を調節する このとき PN Zα < でなく 符号関数 sgn を用いるのは 正しい復元で (a) も十分遠方の点の視線がほぼ平行になり 計算誤差によって に近いところに交点が計 算され Zα となって PN (b) (c) 図 7 図 6 の 3 次元復元結果を上から眺め Google Map の画像に重ねたもの (a) 金谷 三島9) の方法 (b) 本 システムの平均焦点距離の方法 (c) 本システムの固定焦点距離の方法 がそれに引きずられるのを防ぐためである 復元した形状にはスケールの不定性がある これはカメラに近い小さい物体に対してカメ 表 各方法による焦点距離 画素 と再投影誤差 画素 焦点距離のカメラ校正値は 56. 画素 ラを少し移動しても カメラから遠い大きい物体に対してカメラを大きく移動しても 同じ 画像が撮影されるためである6),7) 上記のアルゴリズムでは ktk = となる解を計算して 焦点距離 再投影誤差 いる もしカメラの移動距離 d が既知なら d を あるいはある二つの特徴点の空間の位置 金谷 三島9) 平均焦点距離 固定焦点距離 847.8/85.3.47 33.5.5 48.2.52 例えば r, r 2 の実間隔 d2 が既知なら d2 /kr 2 r k を t とすべての Xα, Yα, Zα に掛け れば実寸の復元となる 第 カメラを基準にするのではなく シーンに固定した X Y Z 世 の抽出と対応付けに SIFT5) を用い RANSA により誤対応を除去した 文献 4) の適 界座標を基準にするには次のように変換すればよい 第 カメラが世界座標の周りに回転 用例を参照 対応付けた特徴点を画像中にマークしている 左はその対応のフロー 画像 行列 R だけ回転し ベクトル t だけ並進した位置にあるとすると 各特徴点のこの世界座 フレームを重ねて対応する特徴点を線分で結んだもの である 図 7(a) はそれらの復元し 標系に関する 3 次元座標 (X α, Y α, Z α ) は次のようになる た 3 次元位置を上方から見たものを Googl Map の画像上に重ねたものである 図 7(b), (c) X α Xα Y α = t + R Yα Z α 7. 実 はそれぞれ本システムの平均焦点距離および固定焦点距離による復元である 表 はカメ (38) ラ校正による焦点距離 OpenV で提供されている方法を用いた およびそれぞれの方 法で得られた値と再投影誤差を示す 単位は画素 Zα 図 6 はほぼ注視画像であり 虚数焦点距離が生じやすいが この例ではどの場合も実数焦 験 点距離が得られている 図 7(a) では自由焦点距離の方法で得られた f, f をそのまま使っ 図 6 の 2 画像から特徴点を抽出し 金谷 三島9) の方法で 3 次元復元を行った 特徴点 ている しかし 値が不正確なために 直交すべき角度がやや鋭角になっている 図 7(b), 7 c 29 Information Processing Society of Japan
(c) 7(b) 7(c) 8. 9) 9) ),2) Taubin 9) 3 9) 2 2 3 6) 7) : (No. 2567, 2572) ) S. Bougnoux, From projective to Euclidean space under any practical situation, a criticisim of self calibration, Proc 6th Int. onf. omput. Vision, January 998, Bombay, India, pp. 79 796. 2) R. I. Hartley, In defense of the eight-point algorithm, IEEE Trans. Patt. Anal. Mach. Intell., 9-6 (997-6), 58 593. 3) R. Hartley and. Silpa-Anan, Reconstruction from two views using approximate calibration, Proc. 5th Asian onf. omput. Vision, January 22, Melbourne, Australia, Vol., pp.338 343. 4) R. I. Hartley and P. Sturm, Triangulation, omput. Vision Image Understand., 68-2 (997-), 46 57. 5) R. Hartley and A. Zisserman, Multiple View Geometry in omputer Vision, ambridge University Press, ambridge, U.K., 2. 6),,, 99. 7) K. Kanatani, Statistical Optimization for Geometric omputation: Theory and Practice Elsevier, Amsterdam, the Netherlands, 996; reprinted, Dover, York, NY, U.S.A., 25. 8),, :,, 2-VIM-2-7 (2-), 49 56. 9),,, :, 42-SIG 6 (2-6), 8. ),, FNS,, 27-VIM-58-4, (27-3), 25 32. ),,,, 27-VIM-6-9, (27-9), 49 56. 2),,, 28-VIM-64- (28-9), 7 24. 3),,, 2 : Hartley vs.,, 28-VIM-62-54 (28-3), 335 342. 4),,, : VIM, 44-SIG 7 (23-2), 7 77. 5) D. G. Lowe Distinctive image features from scale-invariant keypoints, Int. J. omput. Vision, 6-2 (24-), 9. 6),,,,,, 24-VIM-45-2 (24-9), 9 6. 7),,,,, 25-VIM-49-2 (25-5), pp. 7 4. 8),,,, 26-VIM-53-32 (26-3), 27 24. 9),,, 2 3,, 23-VIM-4-2 (23-), 79 86. http://www.img.tutkie.tut.ac.jp/programs.html 8 c 29 Information Processing Society of Japan