【文献】
Jiangsheng You et al.,"Image Matching for translation, rotation and uniform scaling by the Radon Transform",Proceedings 1998 International Conference on Image Processing,米国,IEEE,1998年10月 4日,Vol.1,pp.847-851
【文献】
Yosi Keller et al.,"Pseudopolar-Based Estimation of Large Translations, Rotations, and Scaling in Images",IEEE Transactions on Image Processing,米国,IEEE,2004年12月20日,Vol.14, No.1,pp.12-22
【文献】
Weisi Lin et al,"Perceptual visual quality metrics: A survey",Journal of Visual Communication and Image Representation,NL,Elsevier Inc.,2011年 5月 1日,vol.22, No.4,pp.297-312
【文献】
黒木 秀一、外2名,“フーリエ変換を利用した画像の回転角の測定”,電子情報通信学会論文誌,日本,社団法人電子情報通信学会,1990年 4月25日,Vol.J73-D-II, No.4,pp.590-596
(58)【調査した分野】(Int.Cl.,DB名)
前記ステップ(a)は、前記基準画像Xと前記歪画像Yのうち小さい方をゼロ値の画素で埋込み、前記基準画像Xと前記歪画像Yの水平方向および垂直方向の寸法(mおよびn)を等しくするステップをさらに含む、請求項1に記載の方法。
前記画像評価モジュール処理モジュールは、前記視覚品質評価を、ピーク信号ノイズ比(PSNR)の決定、構造的類似性(SSIM)指標の計算および視覚情報忠実度(VIF)指標の計算のいずれかによって実行するように構成される、請求項14に記載のシステム。
【発明の概要】
【課題を解決するための手段】
【0008】
本発明の目的は、第1の画像Yを、第2の画像Xとアライメントされた(aligned with/そろえられた、又は、位置合わせされた)レジストレーション画像Y
*に処理する画像レジストレーション方法およびシステムを提供することである。
【0009】
本発明の一形態によれば、第1の画像Yを、第2の画像Xとアライメントされたレジストレーション画像Y
*に処理する画像レジストレーション方法であって、前記方法は、
(a)基準画像Xと歪画像Yの間の倍率(スケールファクタ/scale factor)「a」を決定するステップと、
(b)前記歪画像Yに前記倍率「a」でサイズ変更(resize)し、それによって正規化された歪画像Y2を生成するステップと、
(c)前記基準画像Xと前記正規化された歪画像Y2の間の回転角「θ
0」を決定するステップと、
(d)前記正規化された歪画像Y2を回転角「−θ
0」だけ回転し、それによって回転された歪画像Y3を得るステップと、
(e)前記回転された歪画像Y3を、前記基準画像Xと水平方向および垂直方向でアライメントさせるように平行移動することによってレジストレーション画像Y
*を得るステップと、を含む、画像レジストレーション方法が提供される。
上記方法において、前記ステップ(a)は、前記基準画像Xの画素値の和と前記歪画像Yの画素値の和の比を計算することによって前記倍率「a」を決定するステップを含む。
【0010】
前記ステップ(a)は、前記基準画像Xと前記歪画像Yのうち小さいほうをゼロ値の画素で埋込み(padding)、前記基準画像Xと前記歪画像Yの水平方向および垂直方向の寸法(dimension)を等しくするステップをさらに含む。
【0011】
前記ステップ(b)は、前記歪画像Yをセンタリングし、それによってセンタリングされた歪画像Y1を形成し、前記センタリングされた歪画像Y1を前記倍率「a」でサイズ変更し、正規化された歪画像Y2を生成するステップをさらに含む。
【0012】
前記ステップ(c)は、前記基準画像Xおよび前記正規化された歪画像Y2のラドン変換R1、R2をそれぞれ形成し、前記ラドン変換R1およびR2の行の相互相関を計算して前記回転角「θ
0」を決定するステップを含む。
【0013】
前記相互相関を計算するステップは以下のステップを含む。
(i)前記ラドン変換R1の各行と前記ラドン変換R2の各行の間の円周方向相互相関(circular cross−correlation)の組を計算するステップであって、各円周方向相互相関は前記各行の間の回転オフセット角「θ」を規定する、該ステップ。
(ii)各行に対して、最も高い値を有する円周方向相互相関を選択するステップ。
(iii)各行に対して各選択された円周方向相互相関によって規定される回転オフセット「θ」を決定し、前記回転角「θ
0」を、前記決定された回転オフセット「θ」の中央値に等しくなるよう設定するステップ。
【0014】
前記相互相関を計算するステップはさらに以下のステップを含む。
(iv)前記ラドン変換R1と前記ラドン変換R2の1次元(1D)高速フーリエ変換(FFT)を計算するステップ。
(v)前記ラドン変換R1と前記ラドン変換R2の前記1D FFTを組み合わせて中間結果R3とするステップ。
(vi)前記円周方向相互相関を含む第2の中間結果R4を生成するステップであって、前記第2の中間結果R4は前記中間結果R3の要素に基づく、該ステップ。
(vii)中間結果R3の各行において、最も高い値を有する円周方向相互相関を選択するステップ。
(viii)各選択された円周方向相互相関において規定された対応する回転オフセット「θ」を決定し、前記回転角「θ
0」を前記決定された回転オフセット「θ」の中央値に等く設定するステップ。
【0015】
前記組み合わせるステップは、前記ラドン変換R1の前記1D FFTと、前記ラドン変換R2の前記1D FFTの「共役」関数との成分ごとの積を形成するステップを含む。
【0016】
上記の方法において、前記ステップ(vi)はさら、前記第2の中間結果R4において各要素である行に対して、各角「θ」に対応する最大値を決定するステップを含む。前記方法は、前記最大値の中央値を取るステップをさらに含み、前記中央値は回転角「θ
0」に対応する。
【0017】
前記方法は、前処理された画像Y
*の視覚品質評価を実行するステップ、例えば、視覚品質の評価、ピーク信号ノイズ比(PSNR)の決定、構造的類似性(SSIM)指標の計算、視覚情報忠実度(VIF)指標の計算のいずれか一つを実行するステップを有する。
【0018】
本発明の別の形態によれば、第1の画像を処理して、第2の画像にアライメントされたレジストレーション画像Y
*にする画像レジストレーション方法であって、前記方法は、
(a)前記第1の画像と前記第2の画像の間の倍率「a」を決定するステップと、
(b)前記倍率「a」で前記第1の画像をサイズ変更してサイズ変更されたバージョンの前記第1の画像を生成することによって、前記第1の画像と前記第2の画像のサイズを等しくするステップと、
(c)ラドン変換を使用して、サイズ変更したバージョンの前記第1の画像と前記第2の画像の間の回転角「θ
0」を決定するステップと、
(d)サイズ変更したバージョンの前記第1の画像を回転角「−θ
0」だけ回転し、それを前記第2の画像と水平方向および垂直方向でアライメントするように平行移動した後、レジストレーション画像Y
*として出力するステップと、を含む画像レジストレーション方法が提供される。
【0019】
本発明のさらに別の形態によれば、歪画像Yを基準画像Xにアライメントすることによってレジストレーション画像Y
*に処理するシステムであって、前記システムは、
プロセッサと、
前記プロセッサによって実行するためのコンピュータ読み取り可能な命令が記憶されたメモリ装置と、を含み、
倍率「a」を決定する倍率推定モジュールと、
前記倍率「a」を用いて、前記歪画像Yをサイズ変更された歪画像Y2へとサイズ変更する画像サイズ変更モジュールと、
回転角「θ
0」を推定する回転角決定ユニットと、
前記サイズ変更された歪画像Y2を回転角「−θ
0」だけ回転することによって補正された画像Y3を形成する画像回転モジュールと、
前記補正された画像Y3を平行移動して前記レジストレーション画像Y
*を生成する、画像平行移動モジュールと、を形成するシステムが提供される。
【0020】
上記システムにおいて、前記回転角決定ユニットは、
前記サイズ変更された歪画像Y2および前記基準画像Xのラドン変換R1、R2をそれぞれ生成するためのラドン変換モジュールと、
前記ラドン変換R2およびR1から回転角「θ
0」を推定するラドン変換相関モジュールと、をさらに含む。
前記システムは、前記レジストレーション画像Y
*の視覚品質評価を行うための画像評価モジュール処理モジュールをさらに含む。
【0021】
上記システムにおいて、前記画像評価モジュール処理モジュールは視覚品質評価を以下の任意の1つによって行う。
ピーク信号ノイズ比(PSNR)を決定する。
構造的類似性(SSIM)指標を計算する。
視覚情報忠実度(VIF)指標を計算する。
【0022】
前記システムは、
前記歪画像Yをセンタリングしてセンタリングされた画像Y1を生成する画像センタリングモジュールをさらに含み、
前記画像サイズ変更モジュールは、前記倍率「a」を用いた前記センタリングされた画像Y1の前記サイズ変更された歪画像Y2へのサイズ変更を実行するように構成される。
【0023】
前記システムは、
前記歪画像Yからマスキングされた画像Y0を抽出するための画像マスキングモジュールと、前記マスキングされた画像Y0からセンタリングパラメータ「cr」および「cc」を決定するための重心計算モジュールと、をさらに含み、
前記画像センタリングモジュールは、前記センタリングパラメータ「cr」および「cc」を用いて前記歪画像Yから前記センタリングされた画像Y1を生成するように構成される。
【0024】
前記システムは、
前記補正された画像Y3と前記基準画像Xの間のオフセットベクトル「TV」を決定する平行移動推定モジュールをさらに含み、
前記画像平行移動モジュールは、前記補正された画像Y3を前記オフセットベクトル「TV」を用いて平行移動して、前記レジストレーション画像Y
*を生成するように構成される。
【0025】
本発明のさらに別の実施形態のシステムおよび方法は以下に記載される。
【0026】
本発明のさらに別の形態によれば、歪画像Yを基準画像Xにアライメントされたレジストレーション画像Y
*に処理する画像レジストレーション方法であって、
(a)前記基準画像Xと前記歪画像Yとの間の倍率「a」を決定するステップと、
(b)前記基準画像Xを前記倍率「a」の逆数でサイズ変更し、それによって、サイズ変更された基準画像X2を生成するステップと、
(c)前記サイズ変更された基準画像X2と前記歪画像Yの間の回転角「θ
0」を決定するステップと、
(d)前記歪画像Yを回転角「−θ
0」だけ回転することによって、回転された歪画像Y3を決定するステップと、
(e)前記回転された歪画像Y3を前記倍率「a」でサイズ変更し、それによって、前記レジストレーション画像Y
*を得るステップと、を含む画像レジストレーション方法が提供される。
【0027】
上記方法において、前記ステップ(a)は、前記基準画像Xの画素値の和と前記歪画像Yの画素値の和の比を計算することにより倍率「a」を決定するステップを含む。
【0028】
前記ステップ(a)は、前記基準画像Xと前記歪画像Yのうち小さい方をゼロ値の画素で埋込み、前記基準画像Xと前記歪画像Yの水平方向および垂直方向の寸法を等しくするステップをさらに含む。
【0029】
前記ステップ(c)は、前記歪画像Yをセンタリングしてセンタリングされた歪画像Y1を生成するステップと、前記サイズ変更された基準画像X2と前記センタリングされた歪画像Y1の間の回転角「θ
0」の決定を行うステップと、をさらに含む。
【0030】
前記ステップ(c)は、前記サイズ変更された基準画像X2と前記歪画像Yのラドン変換R1およびR2をそれぞれ形成するステップと、前記ラドン変換R1およびR2の行の相互相関を計算して、前記回転角「θ
0」を決定するステップと、を含む。
【0031】
前記ステップ(c)は、前記サイズ変更された基準画像X2と前記センタリングされた歪画像Y1のラドン変換R1およびR2をそれぞれ形成するステップと、前記ラドン変換R1およびR2の行の相互相関を計算して、前記回転角「θ
0」を決定するステップと、を含む。
【0032】
前記相互相関を計算するステップは、
(i)前記ラドン変換R1の各前記行と前記ラドン変換R2の各前記行の間の円周方向相互相関の組を計算するステップであって、各円周方向相互相関は前記各行の間の回転オフセット角「θ」を規定する、該ステップと、
(ii)各行に対して、最も高い値を有する円周方向相互相関を選択するステップと、
(iii)各行に対して各選択された円周方向相互相関によって規定された回転オフセット「θ」を決定し、前記回転角「θ
0」を前記決定された回転オフセット「θ」の中央値に等しく設定するステップと、を含む。
【0033】
前記相互相関を計算するステップは、
(iv)前記ラドン変換R1と前記ラドン変換R2の1次元(1D)高速フーリエ変換(FFT)を計算するステップと、
(v)前記ラドン変換R1と前記ラドン変換R2の前記1D FFTを組み合わせて中間結果R3とするステップと、
(vi)前記円周方向相互相関を含む第2の中間結果R4を生成するステップであって、前記第2の中間結果R4は前記中間結果R3の要素に基づく、該ステップと、
(vii)前記中間結果R3の各行において、最も高い値を有する円周方向相互相関を選択するステップと、
(viii)各選択された円周方向相互相関において規定された対応する回転オフセット「θ」を決定し、前記回転角「θ
0」を前記決定された回転オフセット「θ」の中央値に等しく設定するステップと、をさらに含む。
【0034】
前記組み合わせるステップ(v)は、前記ラドン変換R1の前記1D FFTと、前記ラドン変換R2の前記1D FFTの「共役」関数との成分ごとの積を形成するステップを含む。
【0035】
前記ステップ(vi)は、前記第2の中間結果R4において各行に対して、各角「θ」に対応する最大値を決定するステップを含む。
【0036】
前記方法は、前記最大値の中央値,前記回転角「θ
0」,を選ぶステップをさらに含み、前記回転角「θ
0」は中央値により示された前記相互相関により参照された角度に対応する。
【0037】
前記方法は、前記レジストレーション画像Y
*の視覚品質評価を実行するステップをさらに含む。
【0038】
視覚品質評価が、ピーク信号ノイズ比(PSNR)の決定、構造的類似性(SSIM)指標の計算、視覚情報忠実度(VIF)指標の計算のいずれかにより実行される請求項12の方法。
【0039】
本発明の他の一の観点によれば、第1の画像を第2の画像にアライメントされたレジストレーション画像Y
*に処理する画像レジストレーション方法であって、
(a)前記第1の画像と前記第2の画像の間の倍率「a」を決定するステップと、
(b)前記倍率「a」で前記画像の一方をサイズ変更することによってそれぞれのサイズ変更したバージョンを生成し、前記第1の画像と前記第2の画像のサイズを等しくするステップと、
(c)ラドン変換を使用して、前記第1の画像と前記第2の画像またはそれらの前記対応するサイズ変更したバージョンの回転角「θ
0」を決定するステップと、
(d)前記第1の画像を回転角「−θ
0」だけ回転し、それを前記サイズ変更されたバージョンの前記第2の画像と水平方向および垂直方向でアライメントするように平行移動し、まず、前記倍率「a」でサイズ変更した後、前記レジストレーション画像Y
*として出力するステップと、を含む画像レジストレーション方法が提供される。
【0040】
本発明のさら他の一の観点によれば、歪画像Yを基準画像Xにアライメントすることによってレジストレーション画像Y
*に処理するシステムであって、
プロセッサと、
前記プロセッサによって実行するためのコンピュータ読み取り可能な命令が記憶されたコンピュータ読み取り可能な記憶媒体と、を含み、
前記基準画像Xと前記歪画像Yの間の倍率「a」を決定する倍率推定モジュールと、
前記基準画像Xと前記歪画像Yの間の回転角「θ
0」を推定する回転角決定ユニットと、
前記歪画像Yを回転角「−θ
0」だけ回転することによって回転された歪画像Y3を形成する画像回転モジュールと、
前記回転された歪画像Y3をサイズ変更してレジストレーション画像Y
*を生成する画像スケーリングモジュールを形成するシステムが提供される。
【0041】
前記システムは、センタリングされた画像Y1を生成するための任意的なセンタリングモジュールをさらに含み、前記回転角決定ユニットは、前記基準画像Xと前記センタリングされた画像Y1の間の前記回転角「θ
0」を推定するように構成される。
【0042】
前記システムは、前記倍率「a」の逆数で基準画像Xをスケーリングしてサイズ変更された基準画像X2とする画像前スケーリングモジュールをさらに含み、前記回転角決定ユニットは、前記サイズ変更された基準画像X2と前記センタリングされた画像Y1の間の前記回転角「θ
0」を推定するように構成される。
【0043】
前記システムは、前記補正された画像Y3と前記基準画像Xの間のオフセットベクトル「TV」を決定する平行移動推定モジュールをさらに含み、前記画像平行移動モジュールは、前記回転された歪画像Y3を前記オフセットベクトル「TV」を用いて平行移動させ、レジストレーション画像Y
*を生成するように構成される。
【0044】
前記システムは、レジストレーション画像Y
*の視覚品質評価を実行するための画像評価モジュール処理モジュールをさらに有する。
【0045】
前記画像評価モジュールは、視覚品質評価を以下のいずれかによって実行するように構成される。ピーク信号ノイズ比(PSNR)を決定する。構造的類似性(SSIM)指標を計算する。視覚情報忠実度(VIF)指標を計算する。
【0046】
このようにして、改良された画像レジストレーション方法およびシステムが提供される。
【0047】
本発明の実施形態は、添付の図面を参照して読まれるべき好ましい実施形態の以下の詳細な説明から、より完全に理解される。これらの図面は本願を限定するものとみなされるべきではなく、単に例示を意図したものである。
【発明を実施するための形態】
【0049】
本願発明は、質の測定基準を計算する前に、歪画像を基準画像にアライメントする、歪画像を前処理するシステムおよび方法を提示する。基準画像を歪画像にアライメントする処理(又は、基準画像を歪画像に揃える処理、基準画像を歪画像に位置合わせする処理)は、「レジストレーションすること(registering)」または「レジストレーション(registration)」とも呼ばれることがある。
【0050】
歪画像の画質評価の完全参照(full reference)方法は、歪画像と基準画像との特別な形態の比較に基づく。歪みが2つの画像の間の寸法の違い、その包絡線内での歪画像の水平方向および/または垂直方向の変位または回転などの幾何学的特性の変化、または、他の変形を含む場合、既存の完全参照品質評価方法では、視覚的な評価に対応し得る品質メトリック(測定基準)を生成しない。
【0051】
本発明の実施形態で提示された画像レジストレーション方法の他の応用は、視覚的観察では(ほとんど)同一に見えるが、スケール、回転、または一方の画像にノイズが多いなどの点で異なる、2つの画像の自動比較である。
【0052】
したがって、画像レジストレーションはリモートセンシング(地図作製の更新)、およびコンピュータビジョンに応用される。
【0053】
(異なる時点で撮影された同一の患者の画像データ、例えば、変化の検出や腫瘍の監視など、についての)医療画像レジストレーションが、例えば、呼吸または解剖学上の変化による弾性(非剛性(nonrigid)としても知られている)変形を補償するために用いられ得る。医療画像の非剛性のレジストレーションは、解剖学上のアトラス(atlas)、例えば神経画像のためのタライラッハアトラス(Talairach atlas)などに患者の画像データをレジストレーションするためにも役に立つ。
【0054】
画像レジストレーションのさらなる使用は、天体写真において宇宙を撮影した画像をアライメントすることにある。(自動または手動で入力される)制御点を用いて、一の画像上で主要な特徴を第2の画像にアライメントするように変形を行い得る。
【0055】
また、画像レジストレーションはパノラマ状の画像生成の本質的な部分である。本発明に提示された技術は、カメラおよびカメラつき携帯電話などの組み込み機器においてリアルタイムで実装されて実施されている。
【0056】
図1は、基準画像Xを提供する基準画像源110、歪画像Yを提供する歪み画像源120、および、アライメントパラメータ推定処理140および画像アライメント処理150を含む画像前処理システム130、および、画像評価処理160を有する画像アライメントシステム100の概略ブロック図を示す。
【0057】
基準画像Xおよび歪画像Yの双方が、画像前処理システム130に入力される。ここで、歪画像Yは基準画像Xを参照して前処理され、前処理された画像Y
*を生成し、それはレジストレーション画像Y
*とも呼ばれる。
【0058】
概して、画像評価処理160は入力として、基準画像Xおよび前処理画像Y
*を受信し、評価結果Q0を生成する。
【0059】
歪画像Yは、基準画像Xにさかのぼることができる場合があり、歪画像Yの歪みは、圧縮と、その後の画像の圧縮解除;画像の回転、平行移動、または線形スケーリングによって特徴づけられる幾何学的変形;画像のクロッピングによる画像の一部の喪失;ノイズのあるチャネル上での転送によるノイズの追加のいずれかの単独または組み合わせの結果である場合がある。
【0060】
アライメントパラメータ推定処理140の目的は、基準画像Xを参照して歪画像Yを解析し、歪みの幾何学的特徴の推定を提供するアライメントパラメータAPの組を生成することである。
【0061】
アライメントパラメータAPは歪画像Yを前処理された画像Y
*に変形する画像アライメント処理150によって使用される。前処理された画像Y
*はレジストレーション画像Y
*とも呼ばれることがある。アライメントパラメータAPを用いて、画像アライメントプロセス150は歪画像Yのいくらかまたはすべての幾何学的歪みを補償するが、画像それ自体の視覚的性質は変更しない。結果として、歪画像Yは画像アライメント処理150によってレジストレーション画像Y
*に変形される。すなわち、それは、基準画像Xにレジストレーションされてアライメントされる。
【0062】
画像評価処理160は、従来の視覚品質測定に基づいていてもよく、その場合、評価結果Q0は、構造的類似性(SSIM)指標または視覚情報忠実度(VIF)指標などの品質の測定基準(品質メトリック)を表す。この場合、評価結果Q0は、前処理画像Y
*に関して計算されるが、歪画像Yの真の視覚的品質メトリックとして仮定される。なぜなら、前処理の行為は画像の視覚的内容を変更することなく幾何学的歪みを補正するだけであるためである。例えば、ノイズ及び圧縮の効果は残る。
【0063】
多くの他の応用において、画像評価処理160の出力Q0は、単に、歪画像Yと基準画像Xの間の、例えば両方の画像が同一の場面の別々の写真である場合の、視覚的類似性の指標であり得る。
【0064】
画像前処理システム130の第1の実施形態は、以下に基づく。
−基準画像Xと歪画像Yの両方からキーポイント(keypoint)の組を抽出するステップ。
−最も近い、XおよびYから選択された数のキーポイントのペアをマッチングさせるステップ。
−第1のアフィン変換のパラメータに関して解くステップ。
−第1のアフィン変換の逆(inverse)である第2のアフィン変換を生成するステップ。
−歪画像Yに第2のアフィン変換を適用し、それによって前処理された画像Y
*を生成するステップ。
【0065】
画像内の重要な領域を特徴づけるキーポイントは、D.G.Lowe「Distinctive Image Features from Scale−Invariant Keypoints」、International Journal of Computer Vision、第60巻、第2号、第91−110ページ、2004年によって提示されたスケール不変特徴変換(SIFT)パラメータセットとして知られている。アライメントパラメータ推定処理140の2つの特定の実施形態および画像アライメント処理150の対応する実施形態が提示される。
【0066】
図2aは、画像前処理システム130の第1の実施形態である、アフィン変換に基づく画像アライメントシステム200.1のブロック図を示す。システム200.1は、
図1のアライメントパラメータ推定処理140の実施形態であるアフィン変換パラメータ推定サブシステム210および、
図1の画像アライメント処理150の実施形態であるアフィン画像変換モジュール220を含む。
【0067】
アフィン変換パラメータ推定サブシステム210は、画素埋込モジュール230、X−SIFT抽出モジュール240、Y−SIFT抽出モジュール250、キーポイントマッチングモジュール260、アフィン変換生成モジュール270、および変換反転モジュール280を含む。
【0068】
画素埋込モジュール230においては、画素で計測された水平方向寸法「M」および垂直方向寸法「N」に関して2つの画像XおよびYが同じでないというのでなければ、該2つの画像XおよびYは変更されずに通される。それらの画像が同じサイズでない場合、それらは、画像の一方または両方をゼロ値の画素で適切に埋込むことによって水平方向寸法および垂直方向寸法に関して同じになるように変更される(それぞれの寸法に関して、その寸法(方向)の最小の解像度で画像を埋込む)。寸法を同じにした画像は均等化基準画像X’および均等化歪画像Y’と呼ぶものとする。それらはいずれもM×Nピクセルのサイズである。
【0069】
以下の記載を簡潔にするために、均等化基準画像X’および均等化歪画像Y’は、画素埋込が任意の所与の場合において実際には必要でない場合があるという理解の上で、引き続き、単に基準画像Xおよび歪画像Yと呼ぶことにする。画像処理の概念に詳しい人は、画素埋込が行われた場合、処理のステップ(実施形態の記載には示されない)は、レジストレーション画像Y
*から画素の行または列を削除してレジストレーション画像Y
*の寸法を(真の)基準画像Xの寸法に確実にマッチングさせるように行われ得ることを理解する。
【0070】
X−SIFT抽出モジュール240およびY−SIFT抽出モジュール250において、スケール不変特徴変換(SIFT)の構成は、D.G.Lowe「Distinctive Image Features from Scale−Invariant Keypoints」、International Journal of Computer Vision、第60巻、第2号、第91−110ページ、2004年によって記載されるように、基準画像Xおよび歪画像Yから、それぞれキーポイントの組KxおよびKyの形態で抽出される。
【0071】
キーポイントの組KxおよびKyはキーポイントマッチングモジュール260に入力される。
【0072】
キーポイントマッチングモジュール260において、マッチングしたキーポイントサブセットKx’およびKy’として、所定の数のキーポイントがキーポイントの組「Kx」および「Ky」から選択され、これは、最も近いマッチング特徴ベクトルを有する。選択されたキーポイントサブセットKx’およびKy’は、基準画像Xおよび歪画像Yの間の幾何学的関係を最も近く近似する順アフィン変換「FAT」のパラメータを計算するように、アフィン変換生成モジュール270に送られる。FATは変換反転モジュール280に送られて反転され、逆アフィン変換IATが生成される。これは、
図1のアライメントパラメータの組の第1の実施形態のバージョンを表す。
【0073】
アフィン変換パラメータ推定サブシステム210において、逆アフィン変換IATは歪画像Yをレジストレーション画像Y
*に変換するために適用される。
【0074】
図2bは、画像前処理システム130の第2の実施形態である、代替のアフィン変換に基づく画像アライメントシステム200.2のブロック図である。
図2bの代替のアフィン変換に基づく画像アライメントシステム200.2は、変換反転モジュール280およびアフィン画像変換モジュール220が単一の、より効率的なアフィン補償画像生成モジュール290に結合されている点を除いて、
図2aのアフィン変換に基づく画像アライメントシステム200.1と同じ機能を提供する。
【0075】
図3は、本発明の第1の実施形態によるアフィン変換に基づく再アライメント方法300の流れ図であり、以下のステップを含む。
310 「ゼロ埋込画像」
320 「重要な特徴抽出」
330 「キーポイントをマッチングさせる」
340 「マッチングする画素の数<NK?」
350 「アフィンパラメータを推定する」
360 「変換を反転する」
370 「アフィン変換を実行する」
380 「アフィン補償画像を生成する」、および
390 「空白の画素を埋込む」
【0076】
基準画像Xをx(i,j)(i=1〜m、j=1〜n)として表し、歪画像Yをy(i,j)(i=1〜m’、j=1〜n’)として表す。基準画像Xは、圧縮、ノイズ汚染、ならびにあり得るアフィン変換を含む様々な歪みを受け、歪画像Yとなったものであり得る。2つの画像(XおよびY)が完全に異なる画像であることは考えない。
【0077】
図2aの画像埋込モジュール230において実行されるステップ310「ゼロ埋込画像」において、基準画像Xの寸法mおよびnが、歪画像Yの対応する寸法m’およびn’と比較される。
【0078】
基準画像Xと歪画像Yが同一の数の画像を有するものでない場合、均等化画像XおよびYの寸法をMおよびNに設定する。ここでMはmとm’のうち大きいほうであり、Nはnとn’のうち大きいほうである。いずれかの画像が、オリジナル画像の周囲で寸法の小さいところにゼロ画素を追加することによって埋め込まれて、そのときの場合により、基準画像X及び歪画像Yを置き換える均等化基準画像X’または均等化歪画像Y’を生成し得る。記載を簡潔にするために、埋込が行われたか否かにかかわらず、画像は基準画像Xおよび歪画像Yと呼ぶものとする。
【0079】
図2aのX−SIFT抽出モジュール240およびY−SIFT抽出モジュール250において実施されるステップ320「重要な特徴を抽出する(Extract Key Feature)」において、スケール不変特徴変換(SIFT)は、基準画像Xおよび歪画像Yにおいて、ローカルな(局所的な)特徴をそれぞれ検出および記述するために計算される。詳細な方法は、David Lowe「Distinctive Image Features from Scale−Invariant Keypoints」、International Journal of Computer Vision、第60巻、第2号、第91−110ページ、2004年によって提案されている。
【0080】
SIFTキーポイントが基準画像Xおよび歪画像Yの両方から抽出される。David Loweの「SIFT」手法に続いて、キーポイントのリストが得られ、そこでは各キーポイントはローカルな記述子のベクトル、すなわち画像内での位置、および、色、テクスチャなどのコンテンツの視覚的特徴の記述を含む。基準画像Xおよび歪画像Yに関して、それぞれ複数のKxおよびKyキーポイント(またはその組)が、この方法によって得られる。各キーポイントは128の値からなる記述子ベクトルによって表され、行指標、列指標、倍率、およびキーポイントの方向という4つの値を含む位置レコードと結びつけられてもよい。
【0081】
そして、キーポイント記述子は、N次元ユークリッド空間における点を表すように正規化される。ここでN=128であり、各記述子の128の値に対応する。
【0082】
図2aのキーポイントマッチングモジュール260において実施されるステップ330「キーポイントをマッチングさせる」において、基準画像Xのキーポイントと歪画像Yのキーポイントの間の近接マッチング(close match)の数NKが計算される。ここでマッチングの近さは128次元空間の距離値によって決定される。キーポイントの特徴に関して、その距離が、2番目に近いマッチングのdistRatio倍の距離よりも小さい場合にのみ、マッチングが認められる。パラメータdistRatioは0.2から0.6の範囲内で選択されてよい。distRatioとして0.4の値が以下の実験において用いられる。
【0083】
計算の効率化のために、真のユークリッド距離ではなく単位ベクトル間のドット積を計算するとより効率的であり得る。単位ベクトルのドット積のacosによって計算される角度の比は、小さい角度についてのユークリッド距離の比の良い近似であることに留意する。したがって、最も近い点のベクトル角と2番目に近い近接点のベクトル角の比がdistRatioよりも小さい場合にのみマッチングが保たれる。
【0084】
ステップ340「マッチング画素数<NK?」において、前のステップ330においてマッチングしているキーポイントの対が少なくともNK個見つかったかどうかを判定する。NKの値は好ましくは4に設定される。見つかったマッチング画素の対が4よりも少ない場合(ステップ340からYesで抜ける)、歪画像Yのアフィン前処理を実行することは不可能であり機能は停止しなければならない。そうでなければ(ステップ340からNoで抜ける)、実行はステップ350へと続く。
【0085】
図2aのアフィン変換生成モジュール260で実行されるステップ350「アフィンパラメータを推定する」において、マッチングしているキーポイントのペアの間のグローバルアフィン変換が、線形最小二乗法を用いて決定され得る。
【0086】
詳細には、アフィン変換におけるパラメータは以下のように決定される。
【0087】
基準画像の点[i
kj
k]
Tから歪画像の点[u
kv
k]
Tへのアフィン変換は以下のように記される。
【数1】
ここでkはマッチングしているキーポイントの指数である。未知の変数は変換パラメータa
1、a
2、a
3、a
4、t
1およびt
2である。方程式(1)は以下のように書き直される。
【数2】
ここで、マッチングしているキーポイントの各組が、左端の行列の2つの行と右端のベクトルに寄与する。右端のベクトルはしたがって[u
1v
1u
2v
2...u
NKv
NK]
Tである。この方程式は線形最小二乗法で解くことができる。方程式(2)をPX=Qと略す。したがって、
X=(P
TP)
−1P
TQ (3)
を得る。
【0088】
したがって、変換パラメータはa
1=X(1)、a
2=X(2)、a
3=X(3)、a
4=X(4)、t
1=X(5)、そしてt
2=X(6)である。
【0089】
ステップ350「アフィンパラメータを推定する」の出力は、したがって、上記のようにして導かれた、アフィン変換パラメータa1、a2、a3、a4、t1およびt2の組である。それは基準画像Xと歪画像Yの間の順アフィン変換(FAT)を示す。歪画像は、FATの反転を歪画像に適用することにより、基準画像Xに似せるように、幾何学的に操作することができる。すなわち、スケーリングしたり、回転したり、シフト(平行移動)したり、剪断(shear)したりさえもできる。それによって歪画像Yをレジストレーション画像Y
*、例えばアフィン補償画像y
*(i
*,j
*)に変換し、該画像は基準画像Xと同じ解像度を有する。
【0090】
i
*∈[1,m]かつj
*∈[1,n]として、アフィン補償画像y
*(i
*,j
*)に対して、[i
*j
*]
Tを2つの方法で計算することができる。
【0091】
図2aの変換反転モジュール260において実行されるステップ360「変換を反転する」において、逆アフィン変換は以下のように実施される。
【数3】
その後に、
図2aのアフィン画像変換モジュール220において実行されるステップ370「アフィン変換を実行する」が行われる。それによって前処理された画像が、
y
*(i
*,j
*)=y(u,v) (5)
として得られる。
【0092】
しかしながら、計算された値i
*およびj
*は整数値ではないこともあるため、y
*(i
*,j
*)の要素、整数位置でのアフィン補償画像y
*を見つける問題に取り組む必要がある(存在するが、ステップ370に明示的されていない)。
【0093】
図2bのアフィン補償画像生成モジュール290において実行される、好ましいステップ380「アフィン補償画像を生成する」において、上記のステップ360および370の機能は組み合わせられて、画素が整数指数位置に投影されて直接アフィン補償画像y
*(i
*,j
*)を生成する。画素値はバイリニア(bilinear)補間によって以下のように計算される。
【0094】
アフィン補償画像y
*(i
*,j
*)のすべての画素(i
*,j
*)に対して、以下の式が成立する。
【数4】
【0095】
変形された指標値(u,v)
Tが実数値であるため、バイリニア補間は値y
*(i
*,j
*)を計算するために使用可能である。
【0096】
指標値uおよびvが1≦u≦mかつ1≦v≦nの範囲にある場合、以下の方程式(7)はステップ380の処理の簡略表現を提供する。
【数5】
ここで、
【数6】
および
【数7】
は、それぞれ実数値uの天井関数および床関数である。これに対して、uまたはvがこれらの範囲外である場合、すなわちu<1またはm<uまたはv<1またはn<vである場合、対応する画素値を0に割り当てる:y
*(i
*,j
*)=0。
【0097】
ステップ380、すなわち方程式(6)および(7)の例示的な実施形態は、以下のMATLABコードサンプル1にさらに示される。留意点:MATLAB(matrix laboratory)はMath−Works(http://www.mathworks.com)によって開発された数値計算環境および第4世代プログラミング言語である。MATLABコードサンプル1。
【表1】
【表2】
【0098】
ステップ390「任意的に空白の画素を埋める」において、アフィン補償画像y
*の画像境界近くの任意の空白の画素は空白のままにされてもよく、任意的に、画像境界の非空白画像画素のミラーリング値によって埋込まれてもよい。
【0099】
例えば、以下の画像がある場合、
【表3】
【0100】
この場合、欠けた画素(ゼロ値)のミラーリング値による埋込みは、空白値と非空白値の間の境界における非空白画素値をミラーリングすることによって、実行される。これは以下を与える。
【表4】
【0101】
本発明の第1の実施形態の評価の目的のため、これらの画素は空白のままとする。
【0102】
要するに、アフィン変換に基づく再アライメント方法300は、基準画像Xと歪画像Yを入力として受信し、アフィン変換技術と結合したSIFT技術を用いて、アライメントされたレジストレーション画像Y
*を出力する。
【0103】
図4a、
図4bおよび
図4cは
図3のアフィン変換に基づく再アライメント方法300の画像的な結果を示す。基準画像Xの代わりに用いられるオリジナルの画像が
図4aに示される。アフィン変換に基づく再アライメント方法300を説明するために、
図4aのオリジナル画像のアフィン変換によって、歪画像Yとして使用するための1つの歪画像が作成され、スケーリング、回転、平行移動、および剪断の歪みの組み合わせを呈する。歪画像Yは
図4bに示される。
図4bの歪画像Yは、
図3のアフィン変換に基づく再アライメント方法300により前処理され、前処理画像、すなわち
図4cに示すレジストレーション画像Y
*を生成した。視覚的な検査によって、前処理方法は歪画像の基準画像へのアライメントにおいて優れた仕事を行うことが可能であることがわかる。
【0104】
アフィン変換に基づく再アライメント方法300の性能をテストするために、H.R.Sheikh、Z.Wang、L.CormackおよびA.C.Bovikにより刊行され、http://live.ece.utexas.edu/research/qualityからインターネット上で利用可能な、LIVE画質評価データベースリリース2からの画像を用いて多くの実験が実施された。このデータベースは、5つのタイプの歪みを用いた29のオリジナルの画像に由来する779の歪画像からなる。歪みは、JPEG圧縮、JPEG2000圧縮、ガウスホワイトノイズ(GWN)、ガウスぼかし(GBlur)、およびレイリー高速フェーディング(FF)チャネルモデルを含む。加えて、歪画像は、スケーリング、回転、および水平垂直両方の方向の空間的シフトによってオリジナルの画像から生成された。
【0105】
テストのために、オリジナル(基準)の画像にアライメントされた、前処理されたレジストレーション画像Y
*が、
図1(画像評価プロセス160参照)に示されるように、視覚品質評価をされた。基準画像Xと前処理されたバージョンの歪画像Y、すなわちレジストレーション画像Y
*の視覚的質の測定には、標準的な測定基準が用いられた。3つの完全参照測定基準(full−reference metrics)であるピーク信号ノイズ比(PSNR)、平均構造的類似性指標測定(MSSIM)、および視覚情報忠実度(VIF)測定もまた、画像の視覚的質の測定のために考慮された。
【0106】
3つの性能測定法が考慮された。すなわち、差分平均オピニオン評点(DMOS)と非線形回帰後の客観的モデル出力の間の相関係数(CC)、二乗平均平方根誤差(RMSE)、スピアマン順位相関係数(SROCC)である。SROCCは以下のように定義される。
【数8】
ここでNは評価が行われる画像内の画素の数であり、d
iは主観的評価および客観的評価におけるi番目の画像の順位の差である。Z.WangおよびQ.Li「Information content weighting for perceptual image quality assessment」、IEEE Transactions on Image Processing、第20巻、第5号、第1185−1198ページ、2011年に記載のように、SROCCはノンパラメトリックの順位に基づく相関測定基準であり、主観的評点と客観的評点の間のいかなる単調な非線形マッピングとも独立している。
【0107】
以下の表1−表5は、前処理なし(「オリジナル」方法とする)の歪画像の視覚的質測定と、本発明の第1の実施形態によるアフィン前処理(「提案」方法とする)による歪画像の視覚的質測定との比較から、実験結果をまとめたものである。
【0108】
表1は3つの完全参照測定基準、すなわちPSNR、MSSIMおよびVIFについていかなる歪みも追加されない、このデータベースの実験結果を示す。最良の実施方法は太いフォントで強調する。表1に示された結果は、歪画像の前処理が、オリジナルのデータベースに関するオリジナルの方法よりも画質評価の特性に悪影響を与えないことを示す。すなわち、実際により良くなっている。
【0109】
表2は、LIVEデータベースのすべての画像が、倍率0.9のスケーリング、回転角0.1×180/π度の回転、水平方向に4画素、および垂直方向に4画素の空間シフトにより歪んだ場合の結果を表に示す。この表において、標準の完全参照測定基準、すなわちPSNR、MSSIMおよびVIFが非常に低い質評点を得る。それに対して、前処理方法は非常に良好な視覚的質評点を与える、すなわち、前処理された画像について行われた客観的質評価は、前処理されない画像よりもはるかに主観的評価との相関が高い。最もうまくいく方法は太いフォントで強調する。実際に、前処理なしでは、既存の質測定基準は非常にパフォーマンスが低いことがわかる。
【0110】
表2は、歪みがスケーリング、回転およびシフトを含む場合、オリジナルの方法に対してパフォーマンスが向上することを示す。そのようなパフォーマンスの向上は、
図8、
図9および
図10に示された、提案の測定基準と人の評点(human score)の間の相関がより良いことから見てとれる。
図11、
図12および
図13に示された結果は、本発明の提案の方法はオリジナルの方法よりも、様々な歪みに対してセンシティブではないことを示す。
【0111】
歪画像の提案の前処理に関する目的および根拠が以下に示される。視覚的質にはまったく影響を与えるべきでない、例えば回転、スケーリングなどのアフィン変換を補償することが求められている。事実、
図11および
図12を見ると、提案の方法は小さな回転およびシフトにはかなりフラットである(不変である)ことがわかる。しかし、画像がかなりスケーリングされて、それによって画像の細部が失われる場合、質を同じに保つことが期待できないため、スケーリングはより扱いにくくなる。しかし、VIFに関して
図13cに示すように、この方法は、倍率の変化で緩やかに減少する合理的な評点を与えることが望まれる。この場合には、提案の方法の評点は緩やかに減少して合理的な状態にとどまる(例えばスケーリングが0.5の場合に評点0.7)。これに対してオリジナルの方法(すなわち前処理なし)は、スケーリングが0.5のとき評点が0.1未満となり、オリジナルの解像度の半分に画像をスケーリングした場合に質がそこまで低くはならないであろうため、直感的にあまり意味がない。
【0112】
表1.LIVE画像データベース内のすべての歪テスト画像に関する画質評価の全体的なパフォーマンス比較。
【表5】
【0113】
表2.LIVE画像データベース内の、スケーリング、回転、および水平方向及び垂直方向両方の空間的シフトによって歪みを受けた、各画像に関する画質評価の全体的なパフォーマンス比較。
【表6】
【0114】
表3−表5はPSNR、MSSIMおよびVIFとLIVE画像データベースの5種類の歪みとの対応に関してCC、SROCCおよびRMSEの値をそれぞれ示す。アフィン前変換があるときの評点はそれがないときの評点よりも良い。VIFに関しては、アフィン前処理があるときの評点はそれがないときの評点と同程度である。
表3.LIVE画像データベースに関する非線形回帰後のCC値。
【表7】
表4。LIVE画像データベースに関する非線形回帰後のSROCC値。
【表8】
表5。LIVE画像データベースに関する非線形回帰後のRMSE値。
【表9】
【0115】
図5から
図13はさらに詳細に実験結果を描写する結果の画像を示す。
【0116】
図5a、
図6aおよび
図7aは、アフィン歪みなしの「オリジナル」の測定(すなわち、アフィン変換パラメータ推定およびアフィン補償がない)に関する、LIVE画像データベースのすべての歪画像に対する、DMOS(縦軸)と、それぞれPSNR、MSSIMおよびVIFの評点(横軸)の散布図である。これに対して
図5b、
図6bおよび
図7bは、画像のアフィン前処理後の測定に関する、対応する散布図である。
【0117】
図8a、
図9aおよび
図10aは、追加のアフィン歪み(倍率0.9のスケーリング、回転角0.1×180/π度の回転、水平方向に4画素、および垂直方向に4画素の空間シフト)が加えられた、「オリジナル」の測定に関する、LIVE画像データベースのすべての歪画像に対する、DMOS(縦軸)と、それぞれPSNR、MSSIMおよびVIFの評点(横軸)の散布図である。これに対して
図8b、
図9bおよび
図10bは、画像のアフィン前処理後の測定に関する、対応する散布図である。
【0118】
図8から
図10の散布図は、追加の歪みが加えられた画像の相関を示し、オリジナルの測定と比較しての前処理による顕著な改善を示す。これは水平軸のスケーリングを異なるものにすることにより一層明らかになる。前処理なしの画像の場合、VIF評点は0.004から0.02のような低さである。しかしながら、アフィン前処理後は、VIF評点は0と1の間の範囲にある。
【0119】
図11a、
図11b、および
図11cは、それぞれ測定基準PSNR、MSSIMおよびVIFに関する、LIVEデータベースの選択された画像を平行移動によって歪ませた実験の結果を示すグラフである。平行移動は1から15画素の範囲にある。
【0120】
同様に、
図12a、
図12b、および
図12cは、それぞれ測定基準PSNR、MSSIMおよびVIFに関する、LIVEデータベースの選択された画像を1度から15度までの回転によって歪ませた実験の結果を示すグラフである。
図13a、
図13b、および
図13cは、それぞれ測定基準PSNR、MSSIMおよびVIFに関する、LIVEデータベースの選択された画像を0.5から1.0までの倍率によって歪ませた実験の結果を示すグラフである。
【0121】
図11a、
図11b、
図12a、および
図12bにおける水平軸のスケールは、比較をより明確に示すために(0,0)を含まない。実際、(0,0)においては歪みはなく、したがってオリジナルの画像およびアフィン前処理された画像の両方が同じ測定基準の評点となるであろう。
【0122】
平行移動歪みに関するグラフに見られるように、歪画像をアフィン前処理する提案の方法は、測定基準PSNR、MSSIM、およびVIFに関してほぼ一定の評点を出す。小さな平行移動は視覚的質に影響を与えるべきでないため、ほぼ一定の評点は、良好な前処理が実行される場合に期待される。回転歪みに関しては、提案の方法はMSSIMおよびVIFに関してほぼ一定の評点を出す。PSNR値は空間シフトが増加するにつれて減少するが、提案の方法は前処理なしの標準的なPSNRよりはなお良好である。また、小さな回転は視覚的質に影響を与えるべきでないため、ほぼ一定の評点は、良好な前処理が実行される場合に期待される。測定基準の評点が平行移動および回転においてほとんど一定のままであるという事実は、この測定基準を適用する前の提案の前処理の有効性を示す。スケーリングの歪みに関して、提案の方法はPSNR、MSSIM、およびVIFに関する標準の測定基準よりも顕著に高い評点を出す。倍率が減少するにつれて質の低下が予想されるが、これは、画像の品質は、小さくなるにつれて、低下する(一層ぼやける)ためである。しかしながら、オリジナルの方法において示されるように激しく減少することは直感的にありえない。提案の方法において示されるように、質の緩やかな低下がより適切である。
【0123】
これは
図7a、
図7bおよび
図7cに示されている。そこでは測定基準による値は倍率が減少するにつれて減少する。しかし変化は本発明の前処理方法を用いた場合は緩やかであり、それは、画像のスケーリングに関連した質の緩やかな減少を表す。出願人が発見したことは、画像がアフィン歪みを受ける場合にはMSSIMおよびVIFはPSNRよりも良好な測定基準であるということである。しかしながら、本発明の前処理ステップで生成された3つの測定基準すべてが、前処理なしのこれらの標準的な測定基準よりも高い評点を出す。提案の前処理方法は異なる種類の歪みに対して非常に堅牢であることが容易にわかる。
【0124】
M.P,Sampat、Z.Wang、S.Gupta、A.C.BovikおよびM.K.Markeyらの「Complex wavelet structure similarity:A new image similarity index」、IEEE Transactions on Image Processing、第18巻、第11号、第2385−2401ページ、2009年においてSampat等は、小さい平行移動および小さい回転にのみ不変な、複素ウェーブレット領域における新しい画像類似度指標(CW−SSIM)を提示した。これに対し、本特許出願において提案の前処理方法は、任意量の平行移動、回転、および倍率に対してより堅牢である。このことは、本発明の提案の方法が、異なる種類の歪みのもとでの画像の視覚的な質を測定するのに適していることを示す。
【0125】
画像前処理システム130の第1および第2の実施形態が、SIFTキーポイントマッチングおよびアフィン変換を用いた「アフィン変換に基づく」実施形態であるのに対し、ラドン変換に基づく追加の実施形態が以下に記載される。
【0126】
アフィン変換に基づく方法とラドン変換に基づく方法のいずれも、歪画像を基準画像にアライメントする処理を行うという同じ目的を果たすが、能力にいくつかの相違がある。例えば、アフィン変換に基づく方法はスキューおよびフリップ(反転)を含んだ多くの幾何学的歪みを補正することができるが、ラドン変換に基づく方法はスキューおよびフリップの歪みを補正することができない。その一方で、非常にノイズが多い画像が存在する場合のラドン変換に基づく方法の堅牢性は、アフィン変換に基づく方法の性能を上回る。
【0127】
画像前処理システム130のラドン変換に基づく実施形態は、歪画像Yを基準画像Xによって「レジストレーションする」ための倍率および回転角を決定する方法に基づく。これらの追加実施形態は効率的な実施につながり、ノイズの存在時に優れて堅牢であるが、上記第1および第2の実施形態のほうが優れているある種の歪み(例えば剪断歪み。
図4b参照)に関して制限される。これらの実施形態はさらに、同じ倍率が水平方向と垂直方向の両方にあてはまるという前提に基づいている。しかしながら、提案の方法はノイズが多い環境では非常に良好にはたらく。
【0128】
図14は、以下の機能を実行する複数のブロックを含むラドン変換に基づくアライメントシステムの構成ブロックの組のダイアグラム1400を示す。
【0129】
入力画像、すなわち基準画像Xと歪画像Yの水平方向寸法を均等化する、寸法均等化機能1410。
【0130】
画像をその重心にセンタリングする、センタリング機能1420。
【0131】
歪画像Yのサイズに対する基準画像Xのサイズの比として倍率「a」を決定する、倍率推定機能1430。
【0132】
基準画像Xと歪画像Yの方向の間に存在する回転角θ
0を決定する、角度推定機能1440。
【0133】
歪画像Yまたはそれを倍率「a」によって処理したバージョンをスケーリングして、そのサイズを基準画像Xのサイズにほぼ等しくする、スケーリング機能1450。
【0134】
基準画像Xまたはそれを倍率「a」の逆数によって処理したバージョンをプレスケーリングして、そのサイズを歪画像Yまたはそれを処理したバージョンのサイズにほぼ等しくする、前スケーリング機能1460。
【0135】
歪画像Yまたはそれの処理したバージョンを回転角θ
0だけ回転する、回転機能1470。
【0136】
歪画像Yの処理したバージョンの平行移動を行い、基準画像Xまたはそれの処理したバージョンに垂直方向および水平方向でアライメントする、平行移動機能1480。
【0137】
以下に記載される本発明の実施形態は、歪画像Yと基準画像Xを処理することによってレジストレーション画像Y
*を生成するように、様々な組み合わせおよび順番で構成要素の組1400を使用する。
【0138】
図15は
図1の画像前処理システム130の第3実施形態である、第1のラドン変換に基づく画像アライメントシステム1500のブロック図である。それは以下のモジュールを含む。
【0139】
倍率「a」を決定する、倍率推定モジュール1502。
【0140】
歪画像Yからマスキングされた画像Y
0を抽出する、画像マスキングモジュール1504。
【0141】
マスキングされた画像Y
0からセンタリングパラメータ「c
r」および「c
c」を決定する、重心計算モジュール1506。
【0142】
センタリングパラメータ「c
r」および「c
c」を用いて歪画像Yをセンタリングして、センタリングされた画像Y
1を生成する、画像センタリングモジュール1508。
【0143】
倍率「a」を用いてセンタリングされた画像Y
1をサイズ変更して、サイズ変更された歪画像Y
2とする、画像サイズ変更モジュール1510。
【0144】
サイズ変更された歪画像Y
2は基準画像Xと同じサイズであり、センタリングされているが、さらに回転する必要がある場合がある。
【0145】
第1のラドン変換に基づく画像アライメントシステム1500はさらに以下を含む。
【0146】
サイズ変更された歪画像Y
2と基準画像Xのラドン変換R
1およびR
2をそれぞれ生成するラドン変換モジュール1512、およびラドン変換R2とR1から回転角「θ
0」を推定するラドン変換相関モジュール1516を含む、回転角決定ユニット1511。
【0147】
サイズ変更された歪画像Y
2が回転角「θ
0」だけ回転され、補正された画像Y
3を形成するが、当該画像は基準画像Xからさらに横方向にオフセットしていてもよい、画像回転モジュール1518。
【0148】
補正された画像Y
3と基準画像Xの間にオフセットベクトル「TV」が決定される、平行移動推定モジュール1520。
【0149】
オフセットベクトル「TV」が補正された画像Y
3に適用されてレジストレーション画像Y
*を生成する、画像平行移動モジュール1522。
【0150】
第3の実施形態の方法に関して、歪画像Y内の、実際の視覚的な内容(コンテンツ)の周囲に黒い背景があることが前提となっている。回転角が0度、90度、180度および270度から顕著に離れるにつれ、多くの画素がこの方法においては失われることに留意することも重要である。それらの角度の値から離れるほど、失われる画素の数が多くなる。したがって、アライメント方法の後に画像評価方法が続く場合、これらの角度の値から離れるほど、質の評価はより不正確になる。
【0151】
図16は、第1のラドン変換に基づく画像アライメントシステム1500のモジュールにおいて例示されることができる、第1のラドン変換に基づく方法1600を示し、それは以下のステップを含む。
310 「ゼロ埋込画像」
1610 「倍率aを推定する」
1620 「画像Yの重心を計算する」
1630 「画像YをセンタリングしてY
1を生成する」
1640 「画像Y
1をサイズ変更してY
3を生成する」
1650 「ラドン変換R1およびR2を実行する」
1660 「変換を相関してθ
0を生成する」
1670 「画像Y
1を回転してY
3を生成する」および
1680 「補正された画像Y
3を平行移動する」
【0152】
前記のステップ310「ゼロ埋込画像」で基準画像Xと歪画像Yのサイズを均等化した後、以下の一連のステップが行われ、そこでは、歪画像Yがスケーリングされ、回転され、平行移動され、レジストレーション画像Y
*となる。
【0153】
倍率推定モジュール1502において実行される、ステップ1610「倍率aを推定する」において、フーリエ変換を用いて行うことができる各画像の輝度(luminance)の平均化と、その結果を割ることによって、倍率「a」が決定される。これは、歪画像の生成において、平均の輝度を変更するいかなる操作(例えばコントラストまたはブライトネスの操作など)も行われないことを前提とする。
【0154】
歪画像Yを、基準画像Xを回転してスケーリングしたものとする。
それぞれサイズが(M×N)である、画像XおよびYのフーリエ変換F1およびF2は、以下の方程式(9)及び(10)によって与えられる。
【数9】
【数10】
【0155】
極座標系において、これらの大きさは方程式(11)において与えられるように関連する。B.S.ReddyおよびB.N.Chatterji「An FFT−Based technique for translation,rotation and scale−invariant image registration」、IEEE Transactions on Image Processing、第5巻、第8号、第1266−1271ページ、1996年参照。詳細は、
【数11】
である。
ここでaは倍率(水平方向および垂直方向に同じ倍率を適用する)であり、θ
0は2つの画像の間の回転角であり、
【数12】
θ=tan
−1(n/m) (13)
である。
したがって、倍率は下記のように得られる。
【数13】
【0156】
倍率aは画像XおよびYの2D FFTを計算することなく、より簡潔に得ることもできる。なぜなら
【数14】
および
【数15】
だからである。
【0157】
ここで、歪画像Yがスケーリング同様に平行移動および回転されてもよいと考える。シフト操作が歪画像Yの画素の合計数に対して非常にわずかな数の画素しか取り除かない限り、倍率aを計算するために使用された手法はなおも有効である。
【0158】
画像マスキングモジュール1504において実行されるステップ1620「マスクを作成する」において、マスク画像Y
0が以下の規則によって作成される。Y(m,n)>τの場合にはY
0(m,n)=1、そうでなければY
0(m,n)=0。本発明のこの実施形態において、発明者はτ=40を選んだ。
【0159】
ステップ1630「画像をセンタリングする」において、画像Yはセンタリングされ、センタリングされた画像Y
1となる。マスク画像Y
0の重心(c
r、c
c)は以下のように計算される。
【数16】
【数17】
ここでパラメータc
rおよびc
cは、「重心」とされる画素の、それぞれ行および列の座標を表す。この演算は重心計算モジュール1506で行われる。
【0160】
Y(c
r,c
c)に対応する位置は画像Yの中央に移動され、平行移動は画像センタリングモジュール1508において行われる。センタリングされた画像Y
1は方程式19によって計算される。
【数18】
ここでm∈[1,M]かつn∈[1,N]。(m−c
r+M/2,n−c
c+N/2)が指標範囲[1,M]×[1,N]の外にある場合、Y
1(m,n)は0に設定される。
【0161】
画像サイズ変更モジュール1510において実行される、ステップ1640「画像をスケーリングする」において、倍率「a」は、センタリングされた画像Y
1を正規化された画像Y
2にサイズ変更して、基準画像Xと同じスケールとなるように使用される。これは回転角を推定する準備において行われる。正規化された画像Y
2は、Y
2(m,n)=Y
1(m/a,n/a)と定義される。
【0162】
正規化された画像Y
2の指標mおよびnは整数値でなければならないが、サイズ変更の計算の後は、それらは一般に浮動小数点値であり、切り捨てるか四捨五入して整数値にする必要がある。サイズ変更された画像内の画素値に関して、最も近い画素値を用いてもよく、双一次補間法(bilinear interpolation)を用いて近接の画素の数から画素値を決定してもよい。
【0163】
基準画像Xと正規化された画像Y2の間の回転角θ
0は、以下の2つのステップで対応するラドン変換から得られる。
【0164】
ステップ1650「ラドン変換を実行する」は、ラドン変換モジュール1512において実行され、基準画像Xの対応するラドン変換および正規化された画像Y
2の対応するラドン変換を計算する。
【0165】
基準画像Xのラドン変換R
1は以下のように定義される。
【数19】
そして、正規化された画像Y
2のラドン変換R
2は以下のように定義される。
【数20】
ここでδ()はDiracのデルタ関数である。
【0166】
A.Averbuch、R.R.Coifman、D.L.Donoho、M.Israeli、およびY.Shkolnisky「A framework for discrete integral transformations I − the pseudo−polar Fourier transform」、SIAM Journal on Scientific Computing、30(2)、第764−784ページ、2008において指摘されているように、現代の応用においては離散2Dラドン変換を有することが重要であり、離散2Dラドン変換は最近20年間で多くの著者の注目の対象となっている。最近まで、ラドン変換は2D離散画像のコヒーレントな離散鮮明度を欠いていた、それは代数的に正確で、可逆で、すばやく計算可能である。それゆえ、Averbuch et al.は離散2D画像のための2D離散ラドン変換の観念を定義した。それは絶対値が1未満の勾配の線に沿った合算に基づいている。格子点でない位置における値はゼロ埋込された格子点上の三角補間を用いて定義される。離散ステップがゼロに近づくにつれて、それが連続ラドン変換に収束することから、彼らは、これらの定義が連続体の忠実な記述を提供することを証明した。離散ラドン変換に関するより詳細な記載はAverbuch et al.にある。Averbuch et al.に従い、出願人はそれぞれR
1(r,θ)およびR
2(r,θ)と表される、XおよびY2の離散ラドン変換を行うことができる。それらは共に寸法がK×Lである。
【0167】
2つの画像(基準画像Xおよび正規化された歪画像Y
2)が同一であるか、それらの間の回転およびありうるノイズを除いて少なくとも非常に似ているという前提に基づけば、それらのラドン変換R
1(r,θ)およびR
2(r,θ)は、実質的に同じ、または、θの方向に沿っての円周方向の移動を除いて正確に同じ、すなわちR
2(r,θ)=R
1(r,θ+θ
0)、であるということにつながる。R
1およびR
2のサイズは、半径「r」と角thetaに関して所望の精度に基づいて選択されるパラメータKおよびLによって特徴づけられる。ラドン変換に基づくシステムの実施形態において、パラメータKは用いられるソフトウェア(上記のMATLAB)の「ラドン」関数からデフォルトとして選択される。それはおよそK=sqrt(M*M+N*N)+5、すなわち、画素数において画像の対角線の長さをわずかに上回る画像サイズに関連する。角度の精度に関して、Lの範囲はL=0から179の範囲で選択され、これは約1度の精度を提供する。
【0168】
ラドン変換関数は以下のようにMATLABの文書に記載される。
【0169】
「R=radon(I、theta)は、角度thetaに対して画像強度(intensity image)Iのラドン変換Rを返す。ラドン変換は特定の角度に向けた放射線に沿った画像強度の射影である。thetaがスカラーである場合、Rはtheta度に対するラドン変換を含む列ベクトルである。thetaがベクトルである場合、Rは、各列がtheta内の角度の1つに対するラドン変換である行列である。」
【0170】
理論によれば、ある画像が別の画像を回転したものである場合、該2つの画像は、角度次元に沿った円周方向シフトを除いて、同じラドン変換を有することとなる。今回の場合は、基準画像Xと正規化された歪画像Y
2はより多くの相違点を有する場合があり、変換の同一性は完全には成立しない。この理由のため、いくつかの値rに関して角度を計算し、角度の中央値をとることが、角度の推定をより堅牢なものとする。
【0171】
ラドン変換相関化モジュール1516において実行される、ステップ1660「変換を相関させる」において、2つのラドン変換が相関方法によって体系的に比較され、歪画像Yの基準画像Xに対する最も確からしい傾き角度を推定する。
【0172】
θ
0を計算するために円周方向相互相関を使用し得るが、これは時間を消費する。R
1およびR
2のすべての行に関して、円周方向シフトR
2を行い、その後に該2つの行の相互相関を計算する必要があるため、O(KL
2)の計算複雑性がある。これはO(L
2)の計算複雑性を有する。全部でK個の行があるため、全部の計算複雑性はO(KL
2)である。
【0173】
より効率的な代替例として、1次元(1D)FFTを代わりに用いることもできる。それはO(KLlogL)の計算複雑性を有する。
【0174】
以下の段落はFFTによる相互相関を計算する一般的な方法を提供する。そこには第3の実施形態の好ましい方法がある。
【0175】
出願人は相互相関およびその高速実施の定義をここで簡潔に述べる。2つの離散した実数値関数f[n]およびg[n]の間の相互相関は、D.Nikolic、R.C.Muresan、W.FengおよびW.Singer「Scaled correlation analysis:a better way to compute a cross−correlogram」、European Journal of Neuroscience、第1−21ページ、2012年によれば、
【数21】
である。
f[n]とg[n]の畳み込みはV.Sobolev「Convolution of functions」、Michiel Hazewinkel、Encyclopedia of Mathematics、Springer、ISBN978−1−55608−010−4において以下のように与えられている。
【数22】
関数f[n]およびg[n]の相互相関はh[n]=f[−n]とg[n]の畳み込みと同じである。
【数23】
畳み込み理論に似て、相互相関は以下を満たす。
【数24】
ここでFFTは高速フーリエ変換、conjは複素共役、そしてドットは成分ごとの乗法を意味する。逆FFTをとることによって、f*gの高速実施が可能となる。
【0176】
図17は、ステップ1660「変換を相関させる」の展開を示し、それは以下のステップを含む。
1710「R1およびR2の1D順FFTを計算してR3を生成する」
1720「R3/|R3|に1D逆FFTを適用してR4を生成する」
1730「行に関して最も大きい相互相関係数を決定する」
1740「配列指標最大相互相関係数の中央値からθ
0を推定する」
【0177】
R
1およびR
2のm番目の行をそれぞれr
1(n)およびr
2(n)とし、ここでnは、R1およびR2の円周方向シフトである角度の数、1からLである。
【0178】
ステップ1710「R1およびR2の1D順FFTを計算してR3を生成する」において、1D順FFTはr
1(n)およびr
2(n)のそれぞれに関して計算され、第1の中間結果R3=r
3(u)を生成する。
【数25】
ここで「conj」は複素共役を表し、・は成分ごとの積を表す。
【0179】
ステップ1720「R3/|R3|に1D逆FFTを適用してR4を生成する」において、1D逆FFTはr
3(u)/|r
3(u)|に適用される。すなわち、第2の中間結果R4=r4(n)が計算される。
【数26】
【0180】
r
2(n)=r
1(n+θ
0)であるから、
【数27】
【数28】
r
4(n)=δ(n+θ
0) (24)
したがって、r
4(n)は−θ
0でピーク値1.0をとり、その他のすべての値で0をとる。R4=r4(n)の値のそれぞれがR1とR2の間の、2つのラドン変換の間の与えられたオフセット角に関する相互相関係数のベクトル(すなわち配列)を表す。
【0181】
ステップ1730「行に関して最も大きい相互相関係数を決定する」において、R4=r
4(n)のすべての行から最大値の位置指標θ
*(m)を見つけ出す。すなわち、
【数29】
である。
【0182】
ステップ1740「最大相互相関係数の配列指標の中央値からθ
0を推定する」において、これらの最大値を示すそれぞれの配列指標の中央値をとることにより、2つの画像(基準画像Xとセンタリングされた画像Y
1)の間の回転角の推定をθ
0=median(θ
*(m))として得る。
【0183】
以下のMATLABコードサンプル2が回転角の推定を実施する。
【表10】
【0184】
ここで
図16の記載に戻る。
画像回転モジュール1518において実施されるステップ1670「画像Y
2を回転する」において、サイズ変更された画像Y
2は方向差を補償するように−θ
0度だけ回転され、そして補正された画像Y
3を得る。
【0185】
画像平行移動モジュール1522において実施されるステップ1680「補正された画像Y
3を平行移動する」において、回転された歪画像Y
3は最終位置に平行移動され、それによりレジストレーション画像Y
*を生成する。
【0186】
平行移動のために、G.VargheseおよびZ.Wang「Video Denoising based on a spatiotemporal Gaussian scale mixture model」、IEEE Transactions on Circuits and Systems for Video Technology、第20巻、第7号、第1032−1040ページ、2010年において提示されるグローバルモーション補償(MC)方法を選択してもよい。それは簡潔で、高速で、信頼できる方法であり、整数画素の精度を提供する。Y
3(m,n)は、すでにスケーリングおよび回転に関して補償されている補正された画像とし、X(m,n)を基準画像とする。また、2D FFTをXおよびY
3に実施する。
【数30】
【数31】
そして、相互相関(CC)関数は以下のように定義される。
【数32】
ここでIFFT
2は逆2Dフーリエ変換、・は成分ごとの積、およびconjは複素共役を意味する。推定されるモーションベクトルは以下のように与えられる。
【数33】
平行移動ベクトル推定モジュール1520において計算される、この推定されるモーションベクトルは、回転された歪画像Y
3をレジストレーション画像Y
*としての最終位置に動かすために必要とされる、平行移動ベクトルと解釈される。すなわち、Y
*(m,n)=Y
3(m−m
opt,n−n
opt)。
【0187】
図18−
図22は、以下に記載する
図23において記載された本発明の第4の実施形態において提案されるレジストレーション方法の実行を表す実験の結果を示す。
【0188】
図18a−
図18eから
図22a−
図22eのそれぞれは、
−オリジナルの画像(Original)(a)、
−その歪んだバージョン(Distorted)(b)、
−本発明の第3の実施形態の方法によってレジストレーションされた歪画像(Proposed)(c)、
−SURF方法によってレジストレーションされた歪画像(SURF)(d)、および
−SIFT方法によってレジストレーションされた歪画像(SIFT)(e)
を示す。
【0189】
具体的に、
図18a、
図18b、
図18c、
図18dおよび
図18eはそれぞれ、オリジナルの画像、その歪んだバージョン(スケール=0.8、回転=0.2ラジアン、平行移動=(40,40)、ノイズ標準偏差σ
n=10)、本発明で提案の方法を用いたレジストレーション画像、SURFを用いたレジストレーション画像、SIFTを用いたレジストレーション画像を示す。
【0190】
図19a、
図19b、
図19c、
図19dおよび
図19eはそれぞれ、オリジナルの画像、その歪んだバージョン(スケール=0.8、回転=0.2ラジアン、平行移動=(4,4)、ノイズ標準偏差σ
n=50)、本発明で提案の方法を用いたレジストレーション画像、SURFを用いたレジストレーション画像、SIFTを用いたレジストレーション画像を示す。
【0191】
図18a−
図18eおよび
図19a−
図19eにおいて、3つの方法すべては、比較的小さいものからほどほどのノイズ(それぞれσ
n=10およびσ
n=50)と組み合わされた大きいものからほどほどの空間シフト(それぞれ40,40および4,4)による歪に対して良好に作用することがわかる。
【0192】
図20a、
図20b、
図20c、
図20dおよび
図20eはそれぞれ、オリジナルの画像、その歪んだバージョン(スケール=0.8、回転=0.2ラジアン、平行移動=(4,4)、ノイズ標準偏差σ
n=100)、本発明で提案の方法を用いたレジストレーション画像、SURFを用いたレジストレーション画像、SIFTを用いたレジストレーション画像を示す。この場合には、SURFとSIFTはいずれも失敗であるが、提案の方法は良好に作用する。
【0193】
図21a、
図21b、
図21c、
図21dおよび
図21eはそれぞれ、オリジナルの画像、その歪んだバージョン(スケール=0.5、回転=0.2ラジアン、平行移動=(4,4)、ノイズ標準偏差σ
n=10)、本発明で提案の方法を用いたレジストレーション画像、SURFを用いたレジストレーション画像、SIFTを用いたレジストレーション画像を示す。倍率が0.5では、3つの方法はいずれも良好に作用する。
【0194】
図22a、
図22b、
図22c、
図22dおよび
図22eはそれぞれ、オリジナルの画像、その歪んだバージョン(スケール=0.1、回転=0.2ラジアン、平行移動=(4,4)、ノイズ標準偏差σ
n=10)、提案の方法を用いたレジストレーション画像、SURFを用いたレジストレーション画像、SIFTを用いたレジストレーション画像を示す。倍率が0.1のとき、本発明の実施形態で提案の方法は良好に作用する。しかしながら、SURFおよびSIFTはいずれもうまくいかない。このように、本発明の実施形態で提案の方法は、入力画像が非常に小さくスケーリングされた場合にSIFTおよびSURFよりも良好である。
【0195】
まとめると、σ
n=100のとき、SURFおよびSIFT方法はうまくいかないが、本発明で提案の方法はこの場合に良好に作用する。加えて、倍率a=0.5のとき、3つの方法はすべて良好なレジストレーションを生成する(
図21を参照)。しかしながら、倍率a=0.1のとき、SURFおよびSIFTはいずれもうまくいかないが、本発明で提案の方法は非常に良好に作用する(
図22cを参照)。
【0196】
本発明の第3の実施形態において提案の方法の計算の複雑性が、2つの既存の方法と比較された。従来技術文献において、B.S.ReddyおよびB.N.Chatterji「An FFT−based technique for translation,rotation and scale invariant image registation」、IEEE Transactions on Image Processing、第5巻、第8号、第1266−1271ページ、1996年は、平行移動、回転、およびスケーリングのパラメータを推定する技術を提示する。彼らの方法は6つの順2D FFTおよび6つの逆2D FFTを使用し、時間を消費しエラーが生じやすい。本発明の第3の実施形態では、出願人は平均の比を用いて倍率を推定し、ラドン変換を用いて基準画像と歪画像の間の回転角を推定する。画像のラドン変換は、Y.Pan、K.LiおよびM.Hamdi「An improved constant−time algorithm for computing the radon and hough transforms on a reconfigurable mesh」、IEEE Transactions on Systens,Man, and Cybernetics,Part A−Systems and Humans、第29巻、第4号、第417−421ページ、1999年に示されるように、再構成可能なメッシュ上で一定の時間O(1)で実行可能である。K.Jafari−KhouzaniおよびH.Soltanian−Zadeh「Radon transform orientation estimation for rotation invariant texture analysis」、IEEE Transactions on Pattern Analysis and Machine Intelligence、 第27巻、第6号、第1004−1008ページ、2005年に示されるように、多くとも、ラドン変換は2D FFTの複雑度で実施されることができる。そして出願人は2D FFTを利用してグローバルな空間シフトを得る。G.VargheseとZ.Wang「Video denoising based on a spatiotemporal Gaussian scale mixture model」、IEEE Transacions on Circuits and Systems for Video Technology、第20巻、第7号、第1032−1040ページ、2010年を参照。本発明の第3の実施形態の方法は3つの順2D FFT、2つの順2D FFT、および2つのラドン変換のみを用いる。したがって、出願人の提示する方法は上に引用したB.S.ReddyおよびB.N.Chatterjiの論文に記載の従来技術の方法よりも高速である。E.De CastroおよびC.Morandi「Registration of translated and rotated imaged using finite Fourier transforms」、IEEE Transactions on Pattern Analysis and Machine Intelligence、vol.PAMI−95、第700−703ページ、1987年において、著者は180の逆2D FFTを必要とする画像レジストレーション方法を提示する。そしてさらに、スケーリングの変化が存在する場合はその方法は失敗する。その逆に、本発明の第3の実施形態の方法はスケーリングの変化が存在する場合にも良好に作用し、E.De CastroおよびC.Morandiによって記載される方法よりも計算の複雑性がはるかに低い。
【0197】
本発明の第3の実施形態の方法およびシステムは、既存の従来技術の方法に対して以下の利点を提供する。本発明の第3の実施形態の方法はノイズの多い環境でも非常に良好に作用するが、D.G.Lowe「Distinctive Image Features from Scale−Invariant Keypoints」、International Journal of Computer Vision、第60巻、第2号、第91−110ページ、2004年およびH.Bay、A.Ess、T.TuytelaarsおよびL.Van Gool「SURF:Speeded Up Robust Features」、Computer Vision and Image Understanding (CVIU)、第110巻、第3号、第346−359ページ、2008年によって提示されたものなどの既存の方法は、ノイズのレベルがあまりに高い場合はうまくいかない。加えて、本発明の実施形態の方法は、計算複雑性O(MNlog(MN)+KLlog(L))に関する高速性をもつ。ここで、入力画像はM×Nの寸法であり、ラドン画像はK×Lの寸法である。さらに、本願発明の実施形態の方法は3つのパラメータのすべて(平行移動、回転、およびスケーリング)を推定することが可能である。これに対して大抵の既存の方法は、例えばK.Jafari−KhouzaniおよびH.Soltanian−Zadeh「Radon transform orientation estimation for rotation invariant texture analysis」、IEEE Transactions on Pattern Analysis and Machine Intelligence、 第27巻、第6号、第1004−1008ページ、2005年、または、W.Wei、S.Wang、X.ZhangおよびZ.Tang「Estimation of image rotation angle using interpolation−related spectral signatures with application to blind detection of image forgery」、IEEE Transactions on Information Forensics and Security、第5巻、第3号、第507−517ページ、2010年において提示されたもののように、1つまたは2つのパラメータを計算するのみである。次の章において実施の実験結果は、特にノイズの多い環境において、出願人が提示した方法が画像のレジストレーションに適していることを示す。
【0198】
本発明の第3の実施形態の実験結果
出願人は512×512のBarbara and Lenaの画像においてある実験を行った。そして再びLIVE画質評価データベースリリース2(H.R.Sheikh、Z.Wang、L.CormackおよびA.C.Bovik「LIVE image quality assessment database release 2」、http://live.ece.utexas.edu/research/quality)を使用した。このデータベースは、5つのタイプの歪みを用いて29のオリジナルの画像から作成された779の歪画像を含む。歪みは、JPEG圧縮、JPEG2000圧縮、ガウスホワイトノイズ(GWN)、ガウスぼかし(GBlur)、およびレイリー高速フェーディング(FF)チャネルモデルを含む。出願人はLIVE画像データベース内の779のすべての歪画像において実験を行った。出願人は本発明の第3の実施形態において提示した方法のためにMATLABコードを使用した。
【0199】
表6および表7は、平行移動、回転、およびスケーリング歪みを有し、さらにノイズを加えた、Lena and Barbaraの画像に関する実験的な結果を示す。いずれの表においても、出願人の推定するパラメータは、歪画像を生成した入力歪みパラメータに非常に近い。加えて、出願人の提示する方法はノイズの多い環境において非常に良好に作用する。これに対して、σn=100において、SIFTはマッチングするキーポイントを発見することができず、SURFはLena画像に関して正確な結果を生成しない。
【0200】
本発明の実施形態における提案の方法における回転角の精度が、点Lの数に依存することは、注目すべきである。出願人のシミュレーションにおいて、出願人はラドン変換の回転角においてL=360のサンプル点を用いた。したがって精度は1度(0.0174ラジアンに等しい)である。サンプル点の数を増やせば正確性を上昇させることが可能だが、それにより計算の複雑性が増すことになる。したがって、パラメータの正確性と計算の複雑性の間にはトレードオフの関係がある。
【0201】
表8はLIVEデータベース内のすべての歪画像が、倍率0.9、回転角0.1×180/π度、空間シフト水平方向に4画素、垂直方向に4画素で、さらに歪んだ場合の結果を表にしたものである。このデータベースに関して、ノイズが加えられない場合は、3つの方法(提案のもの、SIFT、およびSURF)すべてがパラメータの良好な推定を可能とする。出願人は、表6および表7に示すように、このデータベース内の画像にノイズを加えた場合、本発明の実施形態の提案の方法がSIFTおよびSURFよりも良好に動作することを推定する。
【0202】
表9は寸法512×512の1対の画像をレジストレーションする3つの方法に関する実行時間を秒で示した表である。本発明の実施形態で出願人が提示する方法はSURFよりも消費時間が短いが、SIFTよりは長いことがわかる。SIFTはC言語で記述されており、したがってMATLABよりも実行が速いことに留意する。H.Bay、A.Ess、T.TuytelaarsおよびL.Van Gool「SURF:Speeded Up Robust Features」、Computer Vision and Image Understanding (CVIU)、第110巻、第3号、第346−359ページ、2008年によれば、仮にSIFTおよびSURFが共にMATLABで記述されていたならば、SIFTはSURFほど速くない。したがって、本発明の実施形態で出願人が提示する方法はSURFおよびSIFTのいずれよりも速いと結論づけることができる。加えて、出願人が実験において示したように、出願人が提示する方法はノイズが多い環境においてSURFおよびSIFTよりも良好に動作する。
【0203】
表6。平行移動された(4,4)、回転された(0.1ラジアン)、およびスケーリングされた(0.9)Lena画像を、異なるノイズレベルでマッチングした結果である。
【表11】
【0204】
表7。平行移動された(4,4)、回転された(0.1ラジアン)、およびスケーリングされた(0.9)Barbara画像を、異なるノイズレベルでマッチングした結果である。
【表12】
【0205】
表8。LIVE画像データベース内のすべてのすでに歪画像に関して、倍率0.9、回転角0.1×180/π度、空間シフト水平方向に4画素、垂直方向に4画素でさらに歪ませたものの、全体として推定されるパラメータを示す。
【表13】
【0206】
表9。サイズ512×512の1対の画像をレジストレーションのための実行時間を示す。SURFおよび提案の方法はMATLABコードを用いるがSIFTはCコードを用いる。H.Bay、A.Ess、T.TuytelaarsおよびL.Van Gool「SURF:Speeded Up Robust Features」、Computer Vision and Image Understanding (CVIU)、第110巻、第3号、第346−359ページ、2008年によれば、仮にSIFTおよびSURFが共にMATLABを用いたならば、SIFTはSURFよりも遅い。
【表14】
【0207】
本発明の第3の実施形態の代替例が、以下に記載される第4の実施形態において開示される。
【0208】
画像前処理システム130の第4の実施形態は第3の実施形態の拡張バージョンであり、歪画像の回転およびサイズ変更の順番が異なっていることを特徴とする。第3の実施形態の動作を要約する。歪画像Yがまずセンタリングされてセンタリングされた画像Y
1となる。そしてサイズ変更されて(スケーリングされて)、基準画像とサイズをマッチングされて、正規化された画像Y
2となる。回転角は基準画像Xおよび正規化された画像Y
2のラドン変換の相関角によって決定される。そして正規化された画像Y
2は回転されて、基準画像Xと方向をアライメントされ、回転された歪画像Y
3となる。そして最後に、回転された歪画像Y
3は水平方向および垂直方向に平行移動されて基準画像Xとアライメントされ、それによってレジストレーション画像Y
*となる。
第4の実施形態は、ラドン変換を計算する前に歪画像Yがスケーリングされるのではなく、第4の実施形態においてはラドン変換を計算する前に基準画像Xが歪画像Yとサイズを同じにされて、歪画像の必要なサイズ変更は回転および平行移動が行われた後で、最後に行われる点で、顕著に異なる。
【0209】
記載の読解を容易にするために、第3の実施形態のいくつかの変数名は第4の実施形態においても再び使用することに留意する。
【0210】
図23は、
図1の画像前処理システム130の第4の実施形態である拡張ラドン変換に基づく画像アライメントシステム2300のブロック図であり、以下のモジュールを含む。
倍率「a」を決定する、
図15の倍率推定モジュール1502。
【0211】
センタリングされた画像Y
1を生成する、オプションのセンタリングモジュール2305。それは
図15の画像マスキングモジュール1504、重心計算モジュール1506、および画像センタリングモジュール1508を含む。
【0212】
倍率「a」の逆数によって基準画像Xをスケーリングして、サイズ変更された基準画像X
2とする、画像前スケーリングモジュール2310。
【0213】
サイズ変更された基準画像X
2およびセンタリングされた画像Y
1それぞれのラドン変換R
1およびR
2を生成するラドン変換モジュール1512、および、ラドン変換R1およびR2から角度を抽出して、角度からシフト効果を消滅させて、回転角「θ
0」を決定する、回転角推定ブロック2320を含む第2の回転角決定ユニット2315。
【0214】
センタリングされた画像Y
1が回転角「θ
0」によって回転されて、回転された歪画像Y
3を形成するが、回転された歪画像Y
3はまだ基準画像Xから横方向にオフセットされている場合のある、
図15の画像回転モジュール1518。
【0215】
回転された歪画像Y
3とサイズ変更された基準画像X
2の間のオフセットベクトル「TV」が決定される、
図15の平行移動推定モジュール1520。
【0216】
オフセットベクトル「TV」が、回転された歪画像Y
3に適用され、補償された画像Y
4が生成される、
図15の画像平行移動モジュール1522。
【0217】
補償された画像Y
4を倍率「a」で基準画像Xのサイズにスケーリングし、そして、レジストレーション画像Y
*を生成する、画像スケーリングモジュール2330。
【0218】
拡張ラドン変換に基づく画像アライメントシステム2300に用いられるモジュール1502、1512、1518、1520および1522は、ここでは、
図15の第1のラドン変換に基づく画像アライメントシステム1500内のものと同じ機能を有するが、それらの入力データおよび出力データは変更されることがある。
【0219】
拡張ラドン変換に基づく画像アライメントシステム2300の、画像Xの前スケーリングの役割を果たすモジュール2310と、画像Y4の最後のスケーリングを行ってレジストレーション画像Y
*を生成する2330は、変更されたトポロジーに適応している。
【0220】
図15の画像マスキングモジュール1504、重心計算モジュール1506、および画像センタリングモジュール1508を含み、同様の機能を有するモジュール2305は、拡張ラドン変換に基づく画像アライメントシステム2300ではオプションである。画像センタリングモジュールの機能は、画像Xを前スケーリングする前に画像Xをセンタリングするために用いられてもよい。
【0221】
拡張ラドン変換に基づく画像アライメントシステム2300の変形において、
図14の構成ブロックは異なって組み合わされてもよい。例えば、回転角を決定する目的で、前スケーリングされた画像Xの代わりにセンタリングされた画像Y
1を前スケーリングすることも可能である。
【0222】
図24は、拡張ラドン変換に基づく画像アライメントシステム2300のモジュールにおいて例示することが可能な、第2のラドン変換に基づく方法2400を示す。それは以下のステップを含む。
【0223】
2410「倍率を推定する」
2420「マスクを生成する」
2430「画像をセンタリングする」
2440「基準画像Xをサイズ変更して、サイズ変更された基準画像X
2を生成する」
2450「回転角を推定する」
2460「歪画像を回転する」
2470「平行移動を推定する」
2480「画像を平行移動する」および
2490「歪画像をサイズ変更する」。
【0224】
基準画像Xおよび歪画像Yが同じ数の画素を有していない場合、前記ステップ310「ゼロ埋込画像」において小さい方の画像の周囲にゼロを埋め込む。ここからは、XおよびYはゼロが埋め込まれたバージョンのXおよびYであるとする。
【0225】
図16のステップ1610と同等のステップ2410「倍率を推定する」において、前記したような倍率「a」が決定される。
【0226】
歪画像Yをセンタリングするためのステップ2420「マスクを生成する」、および2430「画像をセンタリングする」(それぞれ、第1のラドン変換に基づく方法1600の
図16のステップ1620および1630と同等)は、第2のラドン変換に基づく方法2400においてはオプションであり、ここでは必要とされない場合がある。なぜならより堅牢なラドン変換操作が用いられるからである。しかしながら、画像の実際の回転が多くの画素を失わせる(フレームの外側に位置させる)結果となる場合、パラメータの推定は不正確であり得る。したがって、堅牢性を高めるため、重心を計算し、画像をアライメントして(ステップ1620および1630)、歪画像Yをセンタリングされた画像Y
1に変換することを選択してもよい。同様に、ここには示さないが、基準画像Xをオプションとしてセンタリングしてもよい。
【0227】
ステップ2440「基準画像Xをサイズ変更して、サイズ変更された基準画像X
2を生成する」において、計算された倍率「a」の逆数1/aが用いられ、基準画像が歪画像Yと同じスケールとなるように、基準画像Xからサイズ変更された基準画像X
2を生成する。このステップは回転角の推定において有利である。サイズ変更された基準画像をX
2(m,n)=X(m*a,n*a)と表すものとする。
【0228】
ステップ2450「回転角を推定する」において、サイズ変更された基準画像X
2および歪画像Y
1のそれぞれのラドン変換R1およびR2の出力の相関から、基準画像Xおよび歪画像Y
1の間の回転角「θ
0」を推定する。
【0229】
図25はステップ2450「回転角を推定する」の展開であり、以下のステップを有する。
2510「堅牢なラドン変換を計算する」
2520「R1およびR2の1D順FFTを計算してR3を生成する」
2530「1D逆FFTをR3に適用してR4を生成する」
2540「R4の各行において最大値をとる指数(インデックス)を計算する」および
2550「最大値の中央値からθ
0を計算する」。
【0230】
ステップ2510「堅牢なラドン変換R1およびR2を計算する」において、ラドン変換R1およびR2はそれぞれ、サイズ変更された基準画像X
2およびセンタリングされた画像Y
1(ステップ1620および1630における画像のセンタリングがスキップされた場合には、代わりに歪画像Y)に関して計算される。2D離散関数A(x,y)のラドン変換は、A.Averbuch、R.R.Coifman、D.L.Donoho、M.Israeli、Y.Shkolnisky、およびI.SedelnikovのA Framework for Discrete Integral Transformations II−The 2D Discrete Radon Transform、SIAM Journal on Scientific Computing、30(2)、第785−803ページ、2008年に従い、方程式(29)で定義される。
【数34】
ここでδ()はDiracのデルタ関数である。
【0231】
Averbuch et al.によれば、X
2およびY(センタリングがスキップされなかった場合はY
1)の離散ラドン変換を実行可能であり、それぞれR1(r,θ)、R2(r,θ)と表され、両者とも寸法はK×Lである。
【0232】
そしてR1(r,θ)とR2(r,θ)の両者にr方向に沿って1D FFTを行い、得られる行列の絶対値を得ることができる。これは角θの直線に沿うシフトを消去する。例えば、回転の位置が(x0、y0)の場合、角θの直線の空間的シフトは
b=x0cos(θ)+y0cos(θ)
である。
【0233】
この直線に沿ってFFTを行い、大きさを得ることにより、後にシフトbを消去することができる。このことは回転角とグローバルな空間シフトの両者の正確性を向上させる。
【0234】
2つの画像(サイズ変更された基準画像X
2および歪画像YまたはそれのセンタリングされたバージョンY
1)が同一であるか、それらの間の回転およびありえるノイズを除いて少なくとも非常に類似しているという前提に基づくと、それらのラドン変換R
1(r,θ)とR
2(r,θ)は、θ方向に沿う回転シフトを除いて実質的に同一であることがわかる。すなわちR
2(r,θ)=R
1(r,θ+θ
0)である。θ
0を計算するために回転相互相関を用いることができるが、これはO(KL
2)の計算的複雑性で時間を消費する。なぜなら、R
1とR
2のすべての行に関して、R
2を回転シフトして、2つの行の間の相互相関を計算する必要があるからである。これはO(L
2)の計算的複雑性がある。全部でK個の行があるから、全部の計算的複雑性はO(KL
2)となるであろう。
【0235】
より効率的な代替例として、代わりに1次元(1D)FFTを用いることができ、それはO(KLlogL)の計算的複雑性がある。
【0236】
第4の実施形態の好ましい方法が基礎を置くFFTによる相互相関を計算する方法は、第3の実施形態において使用されたものと同一であり、上記ステップ1660の全体的記載を参照すべきであるが、上記シフト効果の消去により増大する。
【0237】
上記考察に基づき、高速相互相関を用いてθ
0を求める。それはO(KLlogL)の計算的複雑性がある。
【0238】
R
1およびR
2のm番目の行をそれぞれr
1(n)およびr
2(n)と表すものとする。ここでn∈[1,L]である。相互相関が最大値になるオフセットが探される。すなわち、2つの行r
1(n)およびr
2(n)が最もマッチングするオフセットが探される。
【0239】
ステップ2520「R1およびR2の1D順FFTを計算してR3を生成する」において、1D順FFTがr
1(n)およびr
2(n)のそれぞれで実行され、第1の中間結果R3=r
3(u)を生成する。
【数35】
ここで「conj」は複素共役を表し、・は成分ごとの積である。
【0240】
ステップ2530「1D逆FFTをR3に適用してR4を生成する」において、1D逆FFTがr
3(u)に適用される。すなわち、第2の中間結果R4=r
4(n)=IFFT(r
3(u))が計算される。第2の中間結果R4は、ラドン変換R1およびR2の行の組み合わせの間の相互相関値のベクトルを構成する。
ステップ2540「R4の各行の中で最大値をとる指数を計算する」において、すべての行に対して、得られた第2の中間結果R4=r
4(n)の最大値の位置指数θ
*(m)を求める。すなわち
【数36】
である。
【0241】
ステップ2550「最大値の中央値からθ
0を計算する」において、θ
0=median(θ
*(m))として、これらの最大値の中から中央値を取る。
以下のMATLABコードサンプル3が回転角の推定を実施する。
【表15】
【0243】
ステップ2460「歪画像を回転する」において、歪画像Yは方向の違いを補償するように−θ度だけ回転される。回転の結果が回転された歪画像Y
3である。
【0244】
ステップ2470「平行移動を推定する」において、サイズ変更された基準画像X
2にアライメントするために必要な、回転された歪画像Y
3の平行移動が、X
2およびX
3の関数として決定される。
【0245】
平行移動のために、G.VargheseとZ.Wang「Video denoising based on a spatiotemporal Gaussian scale mixture model」、IEEE Transacions on Circuits and Systems for Video Technology、第20巻、第7号、第1032−1040ページ、2010年に記載のグローバルモーション補償(MC)方法を選択してもよい。それは簡潔で、高速で、信頼できる方法であり、整数画像の精度を提供する。Y
3(m,n)はすでに回転の補償をされた画像であるとし、X
2(m,m)はスケーリングされた基準画像であるとする。また、
【数37】
【数38】
である。
すると、高速相互相関(FC)関数を方程式(32)において定義できる。
【数39】
ここで、IFFT
2は逆2Dフーリエ変換であり、conjは複素共役であり、記号・は成分ごとの積を示す。推定されたモーションベクトル(平行移動ベクトル)は方程式(33)で与えられる。
【数40】
【0246】
ステップ2480「画像を平行移動する」において、計算された平行移動ベクトルは回転された歪画像Y
3に適用されて平行移動を補償し、平行移動された画像Y
4をもたらす。この関数を形式的にY
4(m,n)=Y
3(m−m
opt,n−n
opt)と示す。
【0247】
ステップ2490「歪画像をサイズ変更する」において、最終的な画像、すなわちレジストレーション画像Y
*を得るために、歪画像Y4が基準画像Xと同じスケールを有するように、計算された倍率aが用いられて歪画像Y4がサイズ変更される。それは今や平行移動、回転、およびスケーリングに関して補償されたものである。この関数を形式的にY
*(m,n)=Y
4(m/a,n/a)と示す。
図23に示す本発明の第4の実施形態のシステムは、部分的にラドン変換のために用いられる処理であるため非常に堅牢であり、第3の実施形態のシステムとは対照的に、いかなる回転角に対しても良好に作用する。
【0248】
図26は、確認の目的で実施される画像アライメントシステム100の実施形態の組み合わせられたブロック
図2600を示す。それは、CPU2604、ネットワークI/Oシステム、およびコマンドインターフェイス2608を含むプロセッサ2602、メモリ、DVD、CD−ROMなどの不揮発性コンピュータ読み取り可能記憶媒体の形態の、前記プロセッサによる実行のためのコンピュータ読み取り可能な命令が記憶されたコンピュータメモリ2610、および、前記画像アライメントシステム100による処理のための基準画像および試験画像を含む画像データベース2612、を含む。前記コンピュータメモリ2610に記憶された前記コンピュータ読み取り可能な命令は、本発明の実施形態の少なくとも1つを実施するためのソフトウェアモジュールを含む。具体的には、アフィン変換に基づく画像アライメントシステム200.1、代替アフィン変換に基づく画像アライメントシステム200.2、第1のラドン変換に基づく画像アライメントシステム1500、拡張ラドン変換に基づく画像アライメントシステム2300である。コンピュータメモリ2610は本発明の実施形態における実施の際に用いられる、画像評価処理160などの画像評価プログラム2614、および、例えばMATLAB環境2618などのソフトウェアライブラリ2614のストレージをさらに有してもよい。
【0249】
図26は、本発明の異なる実施形態の個々のモジュールを、重複するモジュールの組として示している。例えばアフィン変換に基づく画像アライメントシステム200.1と代替アフィン変換に基づく画像アライメントシステム200.2はモジュール230、240、250、260および270を共用する。共用するモジュールに加えて、アフィン変換に基づく画像アライメントシステム200.1はモジュール220および280をさらに含み、代替アフィン変換に基づく画像アライメントシステム200.2はモジュール290を含む。同様に、第1のラドン変換に基づく画像アライメントシステム1500および拡張ラドン変換に基づく画像アライメントシステム2300はモジュールの組を共用する。すなわち1502、1512、1518、1520および1522である。共用するモジュールに加えて、第1のラドン変換に基づく画像アライメントシステム1500はさらにモジュール1504、1506、1508、1510および1516を用い、拡張ラドン変換に基づく画像アライメントシステム2300は2305、2310、2320および2330を含む。
【0250】
本発明の実施形態は、マルチコアCPUであってもよいCPUと、コンピュータ読み取り可能な媒体、例えばメモリ、DVD、CD−ROM、フロッピー(登録商標)、磁気テープ、またはその他の記憶媒体であって、前記CPUによって実行された場合に、上記のようなシステムのモジュールを形成するコンピュータ読み取り可能な命令が記憶されている媒体と、を有する、汎用または専用コンピュータを含む。代替として、
図2a、
図2b、
図15、および
図23は、上記したこれらのシステムのモジュールを形成するために前記CPUによって実行されるコンピュータ読み取り可能な命令を有するコンピュータ読み取り可能な記憶媒体を有する専用コンピュータとファームウェアとの組み合わせまたは特別な専用のハードウェアを含んでもよい。
図2a、
図2b、
図15、および
図23のシステムの各モジュールは、ファームウェアまたは、代替として、
図26に示すようにプロセッサによって実行するための、コンピュータ読み取り可能な記憶媒体内に記憶されたコンピュータ読み取り可能な命令を含んでもよい。
【0251】
結論として、本発明は歪画像を前処理して基準画像にアライメントするための改良されたシステムおよび方法を提供し、それらはより正確な視覚的な画質評価を提供するように使用されてもよく、同様に、歪画像を基準画像にレジストレーションするように作用するように使用されてもよい。
【0252】
本発明の実施形態は詳細に記載されているが、実施形態の変形や変更を以下の特許請求の範囲内で行い得ることは当業者には明らかである。
下記は、本願に当初より記載の発明である。
<請求項1>
歪画像Yを基準画像Xにアライメントされたレジストレーション画像Y
*に処理する画像レジストレーション方法であって、
(a)前記基準画像Xと前記歪画像Yとの間の倍率「a」を決定するステップと、
(b)前記基準画像Xを前記倍率「a」の逆数でサイズ変更し、それによって、サイズ変更された基準画像X2を生成するステップと、
(c)前記サイズ変更された基準画像X2と前記歪画像Yの間の回転角「θ
0」を決定するステップと、
(d)前記歪画像Yを回転角「−θ
0」だけ回転することによって、回転された歪画像Y3を決定するステップと、
(e)前記回転された歪画像Y3を前記倍率「a」でサイズ変更し、それによって、前記レジストレーション画像Y
*を得るステップと、を含む、
画像レジストレーション方法。
<請求項2>
前記ステップ(a)は、前記基準画像Xの画素値の和と前記歪画像Yの画素値の和の比を計算することにより倍率「a」を決定するステップを含む、請求項1に記載の方法。
<請求項3>
前記ステップ(a)は、前記基準画像Xと前記歪画像Yのうち小さい方をゼロ値の画素で埋込み、前記基準画像Xと前記歪画像Yの水平方向および垂直方向の寸法(mおよびn)を等しくするステップをさらに含む、請求項1または2に記載の方法。
<請求項4>
前記ステップ(c)は、前記歪画像Yをセンタリングすることでセンタリングされた歪画像Y1を生成するステップと、前記サイズ変更された基準画像X2と前記センタリングされた歪画像Y1の間の回転角「θ
0」の決定を行うステップと、をさらに含む、請求項1から3のいずれか一項に記載の方法。
<請求項5>
前記ステップ(c)は、前記サイズ変更された基準画像X2および前記歪画像Yのラドン変換R1およびR2をそれぞれ形成するステップと、記回転角「θ
0」を決定するために前記ラドン変換R1およびR2を用いるステップと、を含む、請求項1から4のいずれか一項に記載の方法。
<請求項6>
前記ステップ(c)は、前記サイズ変更された基準画像X2と前記センタリングされた歪画像Y1のラドン変換R1およびR2をそれぞれ形成するステップと、前記ラドン変換R1およびR2の行の相互相関を計算して回転角「θ
0」を決定するステップと、を含む、請求項4または5に記載の方法。
<請求項7>
前記ラドン変換R1およびR2の行の相互相関を計算して回転角「θ
0」を決定するステップは、
(i)前記ラドン変換R1の各行と前記ラドン変換R2の各行の間の円周方向相互相関の組を計算するステップであって、前記各円周方向相互相関は前記各行の間の回転オフセット角「θ」を規定するものである、該ステップと、
(ii)各行に対して、最も高い値を有する円周方向相互相関を選択するステップと、
(iii)各行に対して各選択された円周方向相互相関によって規定される回転オフセット「θ」を決定し、回転角「θ
0」を、前記決定された回転オフセット「θ」の中央値に等しくなるよう設定するステップと、をさらに含む、
請求項5または6に記載の方法。
<請求項8>
前記レジストレーション画像Y
*の視覚品質評価をさらに実行する、請求項1から7のいずれか一項に記載の方法。
<請求項9>
歪画像Yを基準画像Xにアライメントすることによってレジストレーション画像Y
*に処理するシステムであって、
プロセッサと、
前記プロセッサによって実行するためのコンピュータ読み取り可能な命令を有するメモリ装置と、を含み、
前記基準画像Xと前記歪画像Yの間の倍率「a」を決定する倍率推定モジュールと、
前記基準画像Xと前記歪画像Yの間の回転角「θ
0」を推定する回転角決定モジュールと、
前記歪画像Yを回転角「−θ
0」だけ回転することによって回転された歪画像Y3を形成する画像回転モジュールと、
前記回転された歪画像Y3をサイズ変更してレジストレーション画像Y
*を生成する画像スケーリングモジュールと、
を形成するシステム。
<請求項10>
センタリングされた画像Y1を生成するための任意的なセンタリングモジュールをさらに含み、前記回転角決定モジュールは、前記基準画像Xと前記センタリングされた画像Y1の間の回転角「θ
0」を推定するように構成される、請求項9に記載のシステム。
<請求項11>
前記倍率「a」の逆数で前記基準画像Xをスケーリングしてサイズ変更された基準画像X2とする画像前スケーリングモジュールをさらに含み、前記回転角決定モジュールは、前記サイズ変更された基準画像X2と前記センタリングされた画像Y1の間の回転角「θ
0」を推定するように構成される、請求項9または10に記載のシステム。
<請求項12>
前記回転角決定モジュールはさらに、前記サイズ変更された基準画像X2と前記歪画像Yのラドン変換R1およびR2をそれぞれ形成し、回転角「θ
0」を決定するために前記ラドン変換R1およびR2を用いるように構成される、請求項11に記載のシステム。
<請求項13>
前記回転された歪画像Y3と前記基準画像Xの間のオフセットベクトル「TV」を決定する平行移動推定モジュールと、
前記回転された歪画像Y3を前記オフセットベクトル「TV」を用いて平行移動させ、レジストレーション画像Y
*を生成するように構成される画像平行移動モジュールと、
をさらに含む、請求項9から12のいずれか一項に記載のシステム。
<請求項14>
レジストレーション画像Y
*の視覚品質評価を実行するための画像評価モジュール処理モジュールをさらに有する、請求項9から13のいずれか一項に記載のシステム。
<請求項15>
前記画像評価モジュール処理モジュールは、前記視覚品質評価を、ピーク信号ノイズ比(PSNR)の決定、構造的類似性(SSIM)指標の計算および視覚情報忠実度(VIF)指標の計算のいずれかによって実行するように構成される、請求項14に記載のシステム。