第31回信号処理シンポジウム 2016年11月8日 10日 関西大学 位相限定相関法における周波数領域の誤差に着目した 画像の高精度な平行移動量推定 High Accuracy Subpixel Displacement Estimation for Images by Minimizing the

Transcription

1 第3回信号処理シンポジウム 206年月8日 0日 関西大学 位相限定相関法における周波数領域の誤差に着目した 画像の高精度な平行移動量推定 High Accuracy Subpixel Displacement Estimation for Images by Minimizing the Frequency-Domain Error of Phase-Only Correlation 岩田駿人 阿部正英 川又政征 東北大学大学院工学研究科電子工学専攻 Hayato IWATA Masahide ABE Masayuki KAWAMATA Department of Electronic Engineering, Graduate School of Engineering, Tohoku University アブストラクト 本稿では位相限定相関法における周波 かになってきた [4][5] その結果を踏まえると POC 関数 数領域の誤差に着目し 勾配法を用いてこれを最小化す の理論的な形状に基づく既存の手法では 雑音のある環 ることにより 2 つの画像の平行移動量を高精度に推定す 境において必ずしも正確な推定を行うことができないと る手法を提案する 提案法では まず位相限定相関関数 考えられる のピーク座標を検出することでピクセル単位の平行移動 そこで本論文では新たに周波数領域での理論値との誤 量を推定し その値を起点として正規化クロスパワース 差に着目し 勾配法を用いてこれを最小化することによっ ペクトルの理論値との誤差を用いて最急降下法を実行し て 2 つの画像の平行移動量を推定する手法を提案する 提 誤差が最も小さくなるような推定値に収束させる 特に 案法ではまず POC 関数のピークを検出することによりピ 雑音の多い環境において 提案法は従来の手法よりも高 クセル精度での平行移動量を推定し その値を起点とし い精度での推定が可能であること実験により示す て最急降下法により周波数領域での誤差を最小とするよ う値に収束させることで平行移動量を推定する 特に雑 はじめに 音の多い環境において 提案法は従来の手法よりも高い 画像レジストレーション [] とは 2 つの画像の位置合わ 精度での推定が可能であること実験により示す せを行う技術であり 画像処理 映像処理の様々な分野 において基礎となる重要な技術である 2 つの画像の対 応関係を推定するための手法として 近年位相限定相関 法 Phase-Only Correlation, POC[2] に注目が集まって いる POC は高いロバスト性を持つと同時に 後述する相関 ピークへのフィッティングを行うことにより 画像のピク 位相限定相関法を用いた画像の平行移動量推定 2 本章では位相限定相関法を用いて 2 つの画像の平行移 動量を推定する手法について解説する 2. 節で位相限定 相関法の基礎について述べ 2.2 節および 2.3 節 2.4 節 において位相限定相関法を用いてサブピクセル精度での 平行移動量推定を行う既存の各手法について述べる セル間隔よりも細かい精度 サブピクセル精度 での推定 が可能という特徴がある これらの特徴を生かし 動画 2. 位相限定相関法 の動き推定や 3 次元計測など幅広い分野で応用される 位相限定相関法 [2] とは 入力された 2 つの画像に対し 筆者らの研究グループでは 以前より劣化した映像の て周波数領域で振幅の正規化を行った上で相互相関を取 位置ずれを修復するための技術として POC を用いた画像 る手法である 2 つの画像 a および b の周波数スペクト の平行移動量推定に着目してきた 劣化した映像は様々 ルをそれぞれ Am, n, Bm, n としたとき 2 つの画像 なノイズを含むことから高いロバスト性が要求され ま の正規化クロスパワースペクトルは式 で表され これ た最終的に推定値を積算して処理を行うため高い精度で を逆フーリエ変換することにより POC 関数が得られる の推定が必要となる 先行研究において POC 関数を用い ることより高精度な位置ずれ修復が可能であることが示 されている [3] Rm, n = しかし近年 POC 関数の統計的な性質について研究が行 われた結果 雑音のある環境での POC 関数の挙動が明ら Bm, na m, n Bm, nam, n rx, y = F [Rm, n] 2

2 Image a 2.3 二次関数フィッティングを用いた平行移動量推定 PEF を用いた手法より雑音に強い方法として POC 関 数を補間によって拡大し二次関数でフィッティングして ピーク座標を求める手法 POC-QCFM[3] がある この 手法ではピーク近傍の 3 点のみを用いて平行移動量を推 定するため PEF を用いた手法より雑音の影響を受けに くいことが明らかにされている [9] Image b しかし POC-QCFM の推定精度は POC 関数の拡大率 に依存し 平行移動量の値によっては近似による誤差が 図 : POC による画像の平行移動量推定 生じる また 空間領域でのフィッティングに基づいてい ることは変わりないため 雑音によって POC 関数のピー ここで m, n は横方向縦方向それぞれの周波数インデック ス A は A の複素共役 F は逆フーリエ変換をそれ クに変動が生じた場合においては正しく推定できないこ とが予想される ぞれ表している 2 つの画像が循環シフトの関係にある時 POC 関数はシ 2.4 固有値分解を用いた平行移動量推定 上記 2 つの手法とは異なるアプローチで平行移動量を フト量に相当する位置にピークを持つデルタ関数となる 推定する手法として 正規化クロスパワースペクトルか { ら直接平行移動量を推定する手法 [0] がある 具体的に, x = x0, y = y0 3 は特異値分解を用いて正規化クロスパワースペクトルを rx, y = 0, otherwise ランク の行列で近似し 縦横それぞれの特異ベクトル ここで x0, y0 は横方向縦方向それぞれの平行移動量である の位相から最小二乗法によって平行移動量を推定する 一般的に画像の平行移動は循環シフトとは異なるが 図 しかし正規化クロスパワースペクトルから得られる位 に示すように一般的な平行移動の関係においても POC 関 相情報には 2π の任意性があるため unwrap の処理が必 数はデルタ関数状の鋭いピークをもつ よってこのピー 須である 雑音のある環境においてはこの unwrap の処理 ク座標を検出することによりピクセル精度で 2 つの画像 がうまくいかない場合があることから 推定精度が著し の平行移動量を推定することが可能である く低下するという問題がある 2.2 相関ピークの理論式を用いた平行移動量推定 周波数領域の誤差に着目した新しい手法 3 2 つの画像が小数単位の平行移動量を持つとき POC 関数は sinc 関数で近似できるような形状を持つことが知 られている [6] 本章では平行移動量を推定する新しいアプローチとし て周波数領域での理論値との誤差に着目し 勾配法を用い てこれを最小化することにより平行移動量を高精度に推 定する手法を提案する 3. 節において推定値における理 rx, y sincx x0 sincy y0 4 想的な空間シフトと正規化クロスパワースペクトルとの 平均 2 乗誤差 Mean Square Error, MSE を導出し 3.2 この性質を利用し POC 関数の系列に対して sinc 関数を 節において最急降下法を用いてこれを最小化することで フィッティングすることにより x = x0, y = y0 となる本 平行移動量を推定する手法について述べる 来のピーク座標を求め サブピクセル精度で画像の平行移 動量を推定する手法 [7] が知られている また sinc 関数の 性質を利用したピーク評価式 Peak Evaluation Formula, PEF を用いることで ピーク近傍の数点のサンプルから 高い精度で平行移動量を推定する手法 [8] が提案されてお り 高速な処理を可能としている 3. 空間シフトの理論値との誤差 正規化クロスパワースペクトルは 2. 節で述べたよう に式 で表される 2 つの画像が循環シフトの関係にあ るとき それらの正規化クロスパワースペクトルは周波 数領域における空間シフトと等しい しかし PEF は POC 関数が理想的に sinc 関数である 時のみに成立するものである 雑音のある環境において x0 y0 Rm, n = e j2π M m+ N n 5 POC 関数はピークが減衰しサイドローブが上昇するとい う変動をする [5] ため 理想的な相関ピークの形状を前提 ここで x0, y0 は横方向縦方向それぞれの平行移動量 M, N とした手法では正しく推定ができないと考えられる は横方向縦方向それぞれの画像サイズを表している - 8 -

3 正規化クロスパワースペクトル 理想的な空間シフトの系列 Angle 図 2: 正規化クロスパワースペクトルと空間シフト 図 2 に平行移動がある場合の正規化クロスパワースペク トルと理想的な空間シフトそれぞれの偏角を示す 画像端 の影響により位相がわずかに変動しているが Rm, n は おおよそ式 5 に近い形となることが分かる そこでこの Rm, n と理想的な空間シフトの差に着目する Rm, n と推定値 x0, y0 における理想的な空間シフトの MSE εx0, y0 を以下のように定義する εx0, y0 = ここで y0 x0 Rm, n e j2π M m+ N n M N m,n m,n は M /2 m= M /2 N /2 n= N /2 2 6 図 3: Rm, n と理想的な空間シフトとの MSE を表して いる そしてこの εx0, y0 を最小とするような推定値 x 0, y 0 が正しい推定値であると考えられる x 0, y 0 = arg min [εx0, y0 ] M/2, n = N/2 の部分が不連続点となるため この部分を 除外し MSE は以下の式のようになる 7 x0,y0 [ ε x, y = 2 M N { } { } sin M π x sin N π y π y sin π x sin M N 誤差を最小化する手法を考えるために まずこの誤差 の性質を解析する ここでは正規化クロスパワースペク トルが理想的に空間シフトと等しいと仮定し 横方向お よび縦方向における正確な平行移動量と推定値の差をそ 9 れぞれ x = x 0 x0, y = y 0 y0 とする 正規化クロ スパワースペクトルと推定値における理想的な空間シフ トとの MSE は以下のように解くことができる 2 y x e j2π M m+ N n M N m,n [ { }] x y cos 2π m+ n =2 M N m,n M N sinπ x sinπ y 8 = 2 M N sin π x sin π y ε x, y = M N M, N が奇数である場合の ε x, y を図 3 に示す M, N が偶数である場合についても同様の形状となる 図より x, y の範囲におい て ε x, y は単峰性の形状をしており その頂点が x, y = 0, 0 すなわち正しい推定値であることが 分かる 3.2 最急降下法を用いた誤差の最小化 誤差を最小化する問題を考えるとき 一般的には線形 最小 2 乗法などのアルゴリズムが理想的である しかしこ の問題においては 前節において求めたように MSE が非 なお 式 8 は M および N が奇数の場合の式である M 線形の式によって表されるため Rm, n から直接計算に および N が偶数の場合はナイキスト周波数である m = よって平行移動量の推定値 x 0, y 0 を求めることは難しい - 9 -

4 振幅正規化 DFT 振幅正規化 IDFT DFT b 正規化 クロスパワー スペクトル 誤差の計算 a ピーク 検出 a Lenna 最急 降下法 推定値 b Peppers 図 5: 実験に使用した画像 図 4: 最急降下法による画像の平行移動量推定 よって勾配法などの繰り返しの演算によって正しい推 定値を探索する手法を考える 前節において x, y の範囲において MSE が単峰性の形状を 持つことが明らかになったので 正しい推定値が ±[pixel] の範囲に存在することが分かれば勾配法によって収束さ せることが可能である そこで POC 関数のピーク座標検 出を併用することにより平行移動量の推定を行う 図 6: 最急降下法による推定値の収束 図 4 に提案法の概要を示す 提案法ではまず正規化ク ロスパワースペクトルを IDFT することにより POC 関数 を求め そのピーク座標を検出することによりピクセル精 4. 度の平行移動量推定値を求める 次に求めたピクセル精度 の推定値を起点として Rm, n と推定値 x 0 k, y 0 k における理想的な空間シフトとの誤差を用いて最急降下 法を実行し 誤差が最も小さくなるような推定値に収束 させる 提案法の最急降下法における推定値の更新式は 以下の式のようになる x 0 k + y 0 k + = x 0 k y 0 k 提案法における推定値の収束 提案法を用いて画像の平行移動量を推定を行うために まず人工的に位置ずれを発生させた画像を用意した 画像 は図 5 に示す 52 52[pixel] の Lenna を用い フーリ エ変換を用いて横方向に x0 =.3[pixel] 縦方向に y0 = 0.6[pixel] の平行移動を発生させた また 画像端の影響 を低減させるために Hanning 窓による窓かけを行った 初めに 2 つの画像から POC 関数を求めたところ x, y =, の位置にピークが立つ系列が得られた 次 µ εx 0 k, y 0 k に得られたピーク座標を元に推定値の初期値を x 0 0 = 0, y 0 0 = とし 式 0 に基づいてステップサイズ µ = 0. として最急降下法を実行したところ図 6 に示 ここで x 0 k, y 0 k は横方向縦方向の k 回目の更新にお すように正しい推定値 x 0 =.3, y 0 = 0.6 に収束す ける推定値を表しており µ は 回あたりの更新量を調整 ることが確認された 5 回目の更新において推定値は するステップサイズパラメーターである x 0 5 =.2997, y 0 5 = となり 雑音のない 環境においては /00[pixel] を下回る非常に高い精度で の推定が可能であることが分かった 4 実験と比較 本章では 前章で提案した最急降下法を用いる画像の 平行移動量推定の手法について画像を用いた実験を行う 4.2 推定精度の比較 4. 節において実際に最急降下法を用いることで正しい推 定値に収束するかどうかを実験し 4.2 節において画像に 雑音がある環境において既存法および提案法の推定精度 雑音のある環境での提案法および既存の手法について 種類の画像を用いた まず平行移動量あたりの推定精度 を比較する を比較するために x0 = 2 2, y0 = 2 2 の範囲で 推定精度を比較する実験を行った 実験には図 5 に示す 2-0 -

5 PEF based 表 : シミュレーション諸元 フィッティング点数 p=6+6 POC-QCFM Proposed ση2 スペクトル重み幅 U = U2 = 255 拡大率 q=4 更新回数 K = 5 ステップサイズ µ = µ 表 2: 実験に使用したステップサイズ a 画像Lennaを用いた場合 b 画像Peppersを用いた場合 図 7: シフト量あたりの推定誤差 図 8: 雑音分散あたりの推定誤差 人工的に平行移動を発生させた画像を用意し 片方の画 像に分散 ση2 = 0.00 の白色ガウス雑音を付加した なお 化することで 元の画像に付加された雑音の影響を最小 水平方向の平行移動を発生させたときは垂直方向の移動 化できているためだと考えられる 量は y0 = 0 垂直方向の平行移動を発生させた時は水平 次に より一般的な環境を想定し平行移動を x0 = 2 方向の移動量は x0 = 0 とした 生成した二つの画像につ いて 表 に示した諸元に基づき平行移動量の推定を行 2, y0 = 2 2 の範囲でそれぞれ一様分布でランダムに 発生させ 生成した画像のうち片方に ση2 = い 正確な値との誤差を算出した 付加雑音を変えて各 の白色ガウス雑音を付加した シミュレーションをそれ 000 回のシミュレーションを行い 誤差の平均と標準偏 ぞれの雑音分散の値に対して 0000 回行い 推定値の平 差を算出した 画像 Lenna を用いた場合の結果を図 7 に 均 2 乗誤差の平方根 Root Mean Square Error, RMSE 示す 画像 Peppers を用いた場合もほぼ同様の傾向が見 を算出した なおステップサイズ µ は 事前の実験にお られたため結果は省略する いて正しく収束することを確認した表 2 に示すものを用 PEF および QCFM を用いた従来の手法における推定 いた 結果を図 8 に示す 提案法は従来法と比較して雑音を付加した場合の推定 値は平行移動量の値に対して一定の傾向を持っているが 提案法における推定値は正確な値を中心とした分布をし 精度が総合的に向上していることが分かる 特に雑音分 ていることが分かる 提案法は周波数領域の誤差を最小 散 ση2 = 0 2 という厳しい環境においても /00[pixel] - -

6 程度の精度での推定が可能であり 劣化した映像の修復 [6] H. Foroosh, J.B. Zerubia, and M. Berthod, Ex- において効果的に高い精度の推定を行うことが可能であ tension of phase correlation to subpixel registration, IEEE Trans. on Image Processing, vol., no.3, pp , Mar ると考えられる 5 まとめ 本稿では周波数領域における正規化クロスパワースペ クトルと理論値との MSE に着目し 最急降下法を用いて [7] K. Takita, Y. Sasaki, T. Higuchi, and K. Kobayashi, High-accuracy subpixel image registration based に雑音の多い環境下において従来の手法よりも高い精度 on phase-only correlation, IEICE Trans. on Fundamentals of Electronics, Communications and Computer Sciences, vol.86, no.8, pp , Aug. での推定が可能であることを示した これを最小化することでサブピクセル精度で画像の平行 移動量を推定する手法を提案した また実験によって特 今後の検討課題としては まず最急降下法の探索にお [8] S. Nagashima, T. Aoki, T. Higuchi, and K. ける終了条件の検討があげられる 今回の実験では一定 Kobayashi, A subpixel image matching technique using phase-only correlation, 2006 International Symposium on Intelligent Signal Processing and 回数の更新で収束したものと見なしたが 実際の応用に おいては勾配の値などの条件によって探索を終了させる 手法を考える必要がある また 雑音のある環境におい Communications, pp , Dec ては MSE の底となる部分が上昇することから全体として MSE の勾配が緩くなるため 収束に時間がかかってしま うことが分かっている これを解決するために POC 関 数のピークの値などのパラメーター利用してステップサ [9] 松本圭右 八巻俊輔 阿部正英 川又政征 位相限 定相関を用いた平行移動量推定における加法性雑音 イズ µ を適応的に決めてやる必要があると考えられる ソサイエティ大会講演論文集 no.a-4-22 p.89 3 月 の影響の評価 205 年電子情報通信学会基礎 境界 203 参考文献 [] B. Zitova and J. Flusser, Image registration methods: A survey, Image and vision computing, vol.2, no., pp , Oct [0] W.S. Hoge, A subspace identification extension to the phase correlation method, IEEE Trans. on [2] C. Kuglin, The phase correlation image alignment method, Proc. International Conference on Cybernetics and Society, 975, pp.63 65, Jan [3] M. Hagiwara, M. Abe, and M. Kawamata, Estimation method of frame displacement for old films using phase-only correlation, Journal of Signal Processing, vol.8, no.5, pp , Sep [4] S. Yamaki, J. Odagiri, M. Abe, and M. Kawamata, Effects of stochastic phase spectrum differences on phase-only correlation functions part I: Statistically constant phase spectrum differences for frequency indices, 202 3rd IEEE International Conference on Network Infrastructure and Digital Content, pp , Sept [5] 福井一弘 八巻俊輔 阿部正英 川又政征 白色ガウ ス雑音に起因する位相差の変動を持つ実信号の位相 限定相関関数の統計的解析 205 年電子情報通信学 会基礎 境界ソサイエティ大会講演論文集 no.a-4-2 p.6 9 月 Medical Imaging, vol.22, no.2, pp , Feb