3 2 3 2 3 undle Adjustment or 3-D Reconstruction: Implementation and Evaluation Yuuki Iwamoto, Yasuyuki Sugaya 2 and Kenichi Kanatani We describe in detail the algorithm o bundle adjustment or 3-D reconstruction rom multiple images based on our latest research results The main ocus o this paper is the handling o camera rotations and the eiciency o computation and memory space usage when the number o eature points and the number o rames are large An appropriate consideration o these is the core o the implementation o bundle adjustment Doing experiments o undamental matrix computation rom two images and 3-D reconstruction rom multiple images, we evaluate the perormance o bundle adjustment 3 3 3 5),6),2),2) Department o omputer Science, Okayama University 2 Department o omputer Science and Engineering, Toyohashi University o Technolgy 2 3 2 (X, Y, Z) (x, y) X x Y @ y A P () @ Z A P 3 4 (u, v ) t R 4) P = KR / u / I t, K = @ / v / A (2) I K 3) () P X + P 2 Y + P 3 Z + P 4 x = P 3 X + P 32 Y + P 33 Z + P, y = P 2 X + P 22 Y + P 23 Z + P 24 34 (3) P 3 X + P 32 Y + P 33 Z + P 34 P (ij) P ij N (X α, Y α, Z α) M (x α, y α ) ( =,, M, α =,, N) P (3) Vol2-VIM-75 No9 2//2 x, y 3) = 6 c 2 Inormation Processing Society o Japan
E = N M α= = [( pα I α x ) 2 ( α qα + y ) 2 ] α r α r α I α α p α = P X α + P 2 Y α + P 3 Z α + P 4, q α = P 2 X α + P 22 Y α + P 23 Z α + P 24 r α = P 3 X α + P 32 Y α + P 33 Z α + P 34 (5) (x α, y α), α =,, N, =,, M (4) 3 (X α, Y α, Z α ) P 5),6),2),2) 3 3 (X α, Y α, Z α ) P (4) E 3 ( X α, Y α, Z α) P (2), (u, v ) t = (t, t 2, t 3) R,, u, v, t, t 2, t 3 R R 9 R R (4) = I 3 3 R 3 R 3 R R R R + R R = I R R = O ( R R ) = R R R R ω, ω 2, ω 3 ω 3 ω 2 R R = @ ω 3 ω A (6) ω 2 ω 3 ω, ω 2, ω 3 3 SO(3) so(3) 6) a T a T a T 7) (6) ω = (ω, ω 2, ω 3) 6) I ω I (a I)b = a b, (a I)T = a T (6) R R = ω R (7) t t dr /dt ω (7) 7) 7),9) (7) 4 E X α, Y α, Z α, α =,, N,, t, t 2, t 3, u, v, ω, ω 2, ω 3, =,, M 3N + 9M ξ, ξ 2,, ξ 3N+9M ξ k ξ k 2 / ξ k (4) E E N [( )( ) I α pα = 2 xα p α r α r ξ k rα 2 α p α r α ξ k ξ k α= = ( qα + y )( )] α q α r α r α q α r α ξ k ξ k 9) 2 N 2 E ξ k ξ l = 2 α= = ( + I α r 4 α r α q α ξ k [( )( ) p α r α p α r α r α p α r α p α ξ k ξ k ξ l ξ l )] ξ l )( r α q α r α q α r α q α ξ k ξ l ξ k E E/ ξ k 2 2 E/ ξ k ξ l p α, q α, r α p α / ξ k, q α / ξ k, r α / ξ k Vol2-VIM-75 No9 2//2 (8) (9) 2 c 2 Inormation Processing Society o Japan
5 3 (5) p α, q α, r α (X β, Y β, Z β ) δ αβ p α, X β = δ αβ P q α X β = δ αβ P 2, r α X β = δ αβ P 3, p α, Y β = δ αβ P 2 q α 6 Y β = δ αβ P 22, r α Y β = δ αβ P 32, p α = δ αβ P 3 Z β q α Z β = δ αβ P 23 r α Z β = δ αβ P 33 () (2) P P = @ AR I t = @ AK KR I t = @ A u / @ v / AP = u / @ v / AP / = P u P 3 / P 2 u P 32 / P 3 u P 33 / P 4 u P 34 / @ P 2 v P 3 / P 22 v P 32 / P 23 v P 33 / P 24 v P 34 / A () p α, q α, r α λ p α = δ ( λ p α u ) q α r α, = δ ( λ q α v ) r α, λ λ 7 r α λ = (2) (2) P u ( ) ( ) P = @ AR I t = @ AK KR I t u = @ A u / @ v / AP = P 3 P 32 P 33 P 34 @ A (3) / v ( ) P = @ AR I t = @ P v 3 P 32 P 33 P 34 A (4) p α, q α, r α (u λ, v λ ) p α = δ λr α q α r α, =, =, u λ u λ u λ p α q α =, = δ λr α r α, = (5) v λ v λ u λ 8 t (2) P 4 P 4 (R + u R 3 )t + (R 2 + u R 23 )t 2 + (R 3 + u R 33 )t 3 @ P 24 A = KR t = @ (R 2 + v R 3 )t + (R 22 + v R 23 )t 2 + (R 32 + v R 33 )t 3 A P 34 (R 3 t + R 23 t 2 + R 33 t 3 ) (6) P 4 @ P t 24 A = @ t 3 P 34 P 4 @ P 24 P 34 A = @ R + u R 3 R 2 + v R 3 R 3 R 3 + u R 33 R 32 + v R 33 R 33 A, @ t 2 P 4 P 24 P 34 A = @ R 2 + u R 23 R 22 + v R 23 R 23 A A (7) (t λ, t λ2, t λ3 ) tλ (5) tλ p α = δ λ ( r + u r 3 ), tλ p α = δ λ ( r 2 + v r 3 ), tλ p α = δ λ r 3 r, r 2, r 3 r = R R 2 R 3 9, r 2 = R 2 R 22 R 32, r 3 = R 3 R 23 R 33 (2) P Vol2-VIM-75 No9 2//2 (8) (9) 3 c 2 Inormation Processing Society o Japan
P = K(ω R) ω 3 ω 2 ω 2 t 3 ω 3 t 2 I t = KR @ ω 3 ω ω 3 t ω t 3 A (2) ω 2 ω ω t 2 ω 2 t (ω R) = R (ω I) (ω I)t = ω t 7) P / ω, P / ω 2, P / ω 3 R P 3 u R 33 R 2 +u R 23 (t 2 R 3 t 3 R 2 )+u (t 2 R 33 t 3 R 23 ) = @ R ω 32 v R 33 R 22 +v R 23 (t 2 R 32 t 3 R 22 )+v (t 2 R 33 t 3 R 23 ) A, R 33 R 23 (t 2 R 33 t 3 R 23 ) R P 3 +u R 33 R u R 3 (t 3 R t R 3 )+u (t 3 R 3 t R 33 ) = @ R ω 32 +v R 33 R 2 v R 3 (t 3 R 2 t R 32 )+v (t 3 R 3 t R 33 ) A, 2 R 33 R 3 (t 3 R 3 t R 33 ) R P 2 u R 23 R +u R 3 (t R 2 t 2 R )+u (t R 23 t 2 R 3 ) = @ R ω 22 v R 23 R 2 +v R 3 (t R 22 t 2 R 2 )+v (t R 23 t 2 R 3 ) A (2) 3 R 23 R 3 (t R 23 t 2 R 3 ) (ω λ, ω λ2, ω λ3 ) ωλ (5) ωλ p α = δ λ ( r + u r 3 ) (X α t ), ωλ q α = δ λ ( r 2 + v r 3 ) (X α t ), ωλ r α = δ λ r 3 (X α t ) (22) X α = (X α, Y α, Z α ) E 2 E (LM) 9) ( ) X α,, (u, v ), t, R E c = ( 2 ) 2 E/ ξ k, 2 E/ ξ k ξ l, k, l =,, 3N + 9M ( 3 ) ξ k, k =,, 3N + 9M ( + c) 2 E/ ξ 2 2 E/ ξ ξ 2 2 E/ ξ ξ 3N+9M 2 E/ ξ 2 ξ ( + c) 2 E/ ξ2 2 2 E/ ξ 2 ξ 3N+9M @ A 2 E/ ξ 3N+9M ξ 2 E/ ξ 3N+9M ξ 2 ( + c) 2 E/ ξ 3N+9M 2 ξ E/ ξ ξ 2 E/ ξ 2 @ ξ 3N+9M = A @ E/ ξ 3N+9M ( 4 ) X α,, (u, v ), t, R A (23) X α X α + X α, +, (ũ, ṽ ) (u, v ), t t + t, R R(ω )R (24) R(ω ) N [ω ] ω 8) ( 5 ) X α,, (ũ, ṽ ), t, R Ẽ Ẽ > E c c (3) ( 6 ) X α X α,, (u, v ) (ũ, ṽ ), t t, R R (25) Ẽ E δ δ E Ẽ, c c/ (2) (23) c = ξ k 3 R = I, t =, t 22 = (26) 3 2 Y t 2 = 2 Y X Z t 2 = t 23 = (23) ω, ω 2, ω 3, t, t 2, t 3, t 22 3N + 9M 7 (23) LM X α,, (u, v ), t, R (26) (26) X α, t, R X α, t, R X α = ) s R (X α t, R = R R, t = s R (t t ) (27) s = (j, R (t 2 t )) j = (,, ) Vol2-VIM-75 No9 2//2 exp(ω I) so(3) Lie SO(3) 6) 4 c 2 Inormation Processing Society o Japan
2 (8), (9) N MN α= = E/ ξ k (8) ξ k β X β () δ αβ N α= α = β ξ k λ λ, (u λ, v λ ), t λ R λ (2), (5), (8), (22) δ λ = = λ (8) N α= = α α α 2 E/ ξ k ξ l (9) ξ k, ξ l (9) N α= ξ k, ξ l (9) = ξ k, ξ l N α= = E/ ξ k, 2 E/ ξ k ξ l 2 2 E/ ξ k ξ l H(k, l) (3N + 9M 7) (3N + 9M 7) N, M 3N 3 E, 3N 9M F 9M 9 G E 2 E/ X 2 α, 2 E/ X α Y α 2 F 2 E/ X α, 2 E/ X α u 2 G 2 E/ 2, 2 E/ u 2 27NM + 9N + 8M H(k, l) (k, l) 3 3N +9M 7 (23) (3N +9M 7) (3N +9M 7) N, M LU 2) 2 27NM + 6N + 4M I α = 3 (23) (23) E (c) F ( ) ( ) ξ P d P E (c) N F N ξ F F F N G (c) ξ P ξ 3 3N ξ F 9M 7 d P, d F (23) 3N 9M 7 E (c) α, α =,, N α (X α, Y α, Z α) E 2 3 3 (c) ( + c) F α α (X α, Y α, Z α ) α E 2 3 (9M 7) G (c) E 2 (9M 7) (9M 7) ( + c) (28) E (c) E (c) N ξ P + F F N = ξ F = d P, d F ( F (28) ) F N ξ P + G (c) ξ F = d F ξ P 2 ξ F 9M 7 ( N G (c) F α E α (c) F α ) ξ F = α= N α= F α E (c) α α E d F, α E @ (29) E/ X α E/ Y α A (3) E/ Z α ξ F (29) 2 ξ P α X α @ Y α A = E (c) α (F α ξ F + α E) (3) Z α 4 Vol2-VIM-75 No9 2//2 5 c 2 Inormation Processing Society o Japan
2 4) 2745438666 8376652944 768846838 2 76969457 76864356 3 765737 768653 4 778743 76863682 5 76864673 768653 6 7686458 76863682 7 7686458 76866378 σ = LM ɛ LM δ δ = nɛ 2 / 2 n = N α= = Iα ɛ = 2 2 3 2 6 6 = = 6 x, y (Hartley 8 3) ) 22),, 3 2 22) (u, v ) 28 2 2 2 3 9 3 273 (4) e N = 9 E e = (32) N 7 N 7 7 7) σ e 2 / 2 σ 2 N 7 χ 2 N 7 (32) e σ 4) 4) xyx y 3 4 EFNS 3) 2 5 2 4) 22) 3 22) 3 2 3 Oxord http://wwwrobotsoxacuk/~vgg/datahtml 2(a) 36 4983 2 2 P 5266 (23) 2 4 6 2(b) 4983 8 8 3% P (u, v ) t R A 3 A2 n = N Iα 6432 (4) α= = (32) E e = 2n (3N + 9M 7) Vol2-VIM-75 No9 2//2 (33) 6 c 2 Inormation Processing Society o Japan
Vol2-VIM-75 No9 2//2 情報処理学会研究報告 4 ま と め 本論文では多画像からの 3 次元形状復元を行うバンドル調整のアルゴリズムを最新の研究 に基づいて詳細に記述した 本論文で着目したのはカメラ回転の適切な取扱い方 および特 徴点と画像数が多いときの計算とメモリの効率化である まず 2 画像からの基礎行列の計 算に応用し 金谷 菅谷4) の方法と同じ解が得られ 彼らの方法が最適であることを実証 した ただし バンドル調整は効率において彼らの方法に劣る 次に英国 Oxord 大学の実 測データを用いて 3 次元を行った これは特徴点数が非常に多く バンドル調整をその原理 (a) 回 再投影誤差 2 3 4 5 6 7 8 9 327796573463469 2378732275724 7678668765 7232393526 6984294963539 684648452468 67536625569 66882949793228 6642848678532 6639324694876 65756935756945 (b) 回 4 4 42 43 44 45 46 47 48 49 (d) 通りに実装することは困難であるが ここに述べた工夫によって実行できることを示した (c) 理論的にはバンドル調整は任意の 3 次元復元に適用できる万能手法であるが 椋木ら9) 再投影誤差 626388763577 626973343624 62679434579 6264995753774 626262285242 6259944542568 6259624742569 62593353669423 6259486639 6258762887785 が検討したような単純な形状 例えば球面 と単純なカメラ配置 例えば同心円上 から出 発して複雑な形状に収束させることには無理があり 初めに精度のよい近似的な 3 次元復 元を行うことが必要である 近似的 3 次元復元の代表はアフィンカメラ近似による因子分 解法2),5),),2) であるが 透視投影モデルによる自己校正法),),8) のほうが精度が高い バンドル調整はそのような高精度の復元を行った後の最終的な補正を行う手段とみなすのが よいと思われる 謝辞: バンドル調整についてご助言頂いた東北大学の岡谷貴之准教授に感謝します 本研究の一部は文 部科学省科学研究費基盤研究 ( 2572) の助成によった 参 (e) 考 文 献 ) ハノ アッカーマン, 新妻弘崇, 金谷健一, 自己校正法のための射影復元の計算量削減, 情報処理学会研究報告, 27-VIM-6- (27-9), 63 7 2) 浅原 清太郎, 金谷 健一, 菅谷 保之, ハノ アッカーマン, 未校正因子分解法による 3 次元復元 比較実験, 情報処理学会研究報告, 25-VIM-5-2 (25-), 45 52 3) R I Hartley, In deense o the eight-point algorithm, IEEE Trans Patt Anal Mach Intell, 9-6 (997-6), 58 593 4) R Hartley and A Zisserman, Multiple View Geometry in omputer Vision, 2nd ed, ambridge University Press, ambridge, UK, 24 5) 金出武雄, コンラッド ポールマン, 森田俊彦, 因子分解法による物体形状とカメラ運 動の復元 電子情報通信学会論文誌 D-II, J74-D-II-8 (993-8), 497 55 6) K Kanatani, Group-Theoretical Methods in Image Understanding, Springer, erlin, Germany, 99 7) K Kanatani, Statistical Optimization or Geometric omputation: Theory and Practice Elsevier, Amsterdam, the Netherlands, 996; reprinted, Dover, York, NY, USA, 25 8) 金谷健一, 形状 AD と図形の数学, 共立出版, 998 図 2 (a) 用いた 36 画像中の一つ (b) 点を抜き出した場合のヘッセ係数行列の非零パタン (c) 反復回数に対する再投影誤差 e のグラフ (d) 反復回数に対する再投影誤差 e の値 (e) 3 次元 復元 赤 初期位置 緑 収束位置 最小二乗法による初期復元の再投影誤差は e = 327797 画素であったが バンドル調整の 結果 e = 625876 画素に低下した 反復回数による e の変化を図 2(c) に その具体的な数 値を図 2(d) に示す 反復回数は 49 回であり 実行時間は 2 分 5 秒であった ただし ++言語を用い PU には Intel ore2duo E675, 266GHz 主メモリ 4G OS には Windows Vista を用いた 図 2(e) は復元した各点の 3 次元位置をある方向から眺めたもの である 赤は最小二乗法による初期位置であり 緑が最終的な復元位置である 7 2 Inormation Processing Society o Japan
9),,, 25 ) ( ),, 2,, 29, pp 62-68 ),,,,,, 25-VIM-5-6 (25-9), 3 38 2),,,, PRMU23-8 (23-), 9 24 3),, FNS,, 27-VIM-58-4, (27-3), 25 32 4) K Kanatani and Y Sugaya, ompact undamental matrix computation, IPSJ Trans omput Vis Appl 2 (2-3), 59 7 5) M I A Lourakis and A A Argyros, Is Levenberg-Marquardt the most eicient optimization algorithm or implementing bundle adjustment?, Proc th Int on omput Vis, Vol 2, October 25, eijing, hina, pp 526 53 6) M I A Lourakis and A A Argyros, SA: A sotware package or generic sparse bundle adjustment, AM Trans Math Sotware, 36- (29-3), 2: 3 7),,, 3,, 44- (23-), 2864 2872 8),,, :,, 27-VIM-57-5 (27-), 9 6 9),,,,, 4-SIG3 (24-2), 64 73 2),,, 29-VIM-67-37 (29-6), 6 2) Triggs, P F McLauchlan, R I Hartley, and A Fitzgibbon, undle adjustment A modern synthesis, in Triggs, A Zisserman, and R Szeliski, (eds), Vision Algorithms: Theory and Practice, Springer, erlin, 2, pp 298 375 22),,,, 3,, 29-VIM-68-5 (29-9), 8 A ( P = Q ) q P 3 3 Q 4 q () P det Q < Q q Q = ckr, q = ckr t (34) c t t = Q q (35) R R R = I (34) QQ = c 2 KR RK = c 2 KK (36) (QQ ) = c 2 (K ) (K ) (37) (QQ ) (QQ ) = (38) (37), (38) = c K = ck (39) (34) (39) Q = R (4) R R = (Q) (4) K (39) (3,3) A2 3 (3) xp 3 X + xp 32 Y + xp 33 Z + xp 34 = P X + P 2 Y + P 3 Z + P 4 yp 3 X + yp 32 Y + yp 33 Z + yp 34 = P 2 X + P 22 Y + P 23 Z + P 24 (42) p α n α (= I = α) p α 3 X α 2n α x α P 3 P x α P 32 P 2 x α P 33 y α P 3 P 2 y α P 32 P 22 y α P 33 P 23 P 3 X α Y α Z α x α P 34 P 4 = y α P 34 P 24 9) Vol2-VIM-75 No9 2//2 (43) 8 c 2 Inormation Processing Society o Japan