(58)【調査した分野】(Int.Cl.,DB名)
前記状態遷移を、時刻(n−1)におけるシステム状態量と時刻nにおけるシステム状態量との相関を表す忘却係数を用いて表した請求項1または請求項2記載の画像復元装置。
前記状態遷移を、時刻(n−1)におけるシステム状態量と時刻nにおけるシステム状態量との相関を表す忘却係数を用いて表し、前記忘却係数は、時刻(n−1)におけるシステム状態量を構成する表面画像のフレームの要素より裏面画像のフレームの要素に対する値を小さくした請求項4記載の画像復元装置。
【発明を実施するための形態】
【0029】
以下、図面を参照して本発明の実施の形態を詳細に説明する。
【0030】
<第1の実施の形態>
まず、第1の実施の形態について説明する。第1の実施の形態では、動画像を対象として、各フレームの画像を復元する画像復元装置に本発明を適用した場合について説明する。
【0031】
図1に示すように、第1の実施の形態に係る画像復元装置10は、CPU、ROM、RAM、及びHDDを含むコンピュータで構成されており、入力インタフェース部20、操作部30、記憶部40、画像復元処理部50、及び出力インタフェース部60を備えている。
【0032】
入力インタフェース部20は、画像入力装置から提供される画像データをコンピュータで処理可能なデータ形式に変換するなどの入力処理を行う。
【0033】
画像入力装置は、復元処理の対象となる画像データ(劣化画像)をデジタルデータとしてコンピュータに入力するための入力装置であり、動画像を構成する時系列に連続した各フレームの各画像データを画像復元装置10に入力する。画像入力装置としては、例えば、カメラ、記録メディア、モデムなどを用いることができる。カメラは、撮像機能を有する全ての装置を意味し、例えば、デジタルビデオカメラの他に、カメラ機能を搭載した携帯電話や、防犯カメラ(監視カメラ)、画像診断を行うための医療機器(内視鏡、レントゲン、エコー、CT、MRIなど)などを含むことができる。記録メディアは、画像データを記録可能な記録媒体を広く意味し、例えば、磁気ディスク(HDDやFDなど)や、光ディスク(CDやDVD、BDなど)、光磁気ディスク(MO)、フラッシュメモリ(メモリカードやUSBメモリなど)などである。モデムは、外部の通信ネットワーク(例えば、電話回線やLAN、インターネットなど)と接続するための装置である。なお、画像入力装置の種類は、画像復元装置10の用途に応じて適用可能なものを選択すればよい。
【0034】
図示は省略するが、入力インタフェース部20は、画像入力装置の種類に応じて別個独立に設けられる。例えば、記録メディアの入力インタフェース部20は、ドライブとも呼ばれ、記録メディアの種類に応じて様々な種類のドライブが使用可能である。なお、ドライブは、記録メディアを読み書きする装置であり、記録メディアに関する限り、通常は、入力インタフェース部20と出力インタフェース部60とは一体化されている。また、通常、モデムは、画像入力装置としても画像出力装置としても機能しうるため、モデムに関しても、通常、入力インタフェース部20と出力インタフェース部60とは一体化されている。入力インタフェース部20は、コンピュータ本体の内部に格納された内蔵カード(ボード)であってもよいし、内部インタフェース部(図示省略)を介して接続された外部設置型機器であってもよい。
【0035】
なお、画像入力装置が画像情報をアナログデータとして出力する場合、対応する入力インタフェース部20は、サンプリング部及びA/D変換部(共に図示せず)を有する。サンプリング部は、所定のサンプリング周波数で、入力されたアナログ信号をサンプリング処理し、A/D変換部に出力する。サンプリング周波数は、復元処理対象(情報源)の種類に応じて変更可能である。A/D変換部は、サンプリングされた信号の振幅値を所定の分解能でA/D変換処理する。
【0036】
操作部30は、例えば、キーボートやマウス、タッチパネル等であるが、音声認識装置などを用いてもよい。使用者は、操作部30を用い、本装置を操作することができる。また、操作部30は、パラメータ設定部31、領域指定部32、及び処理フレーム指定部33を有する。パラメータ設定部31は、使用者の入力操作により、本実施の形態における画像復元処理に必要な各種パラメータの値を設定する。領域指定部32は、使用者の入力操作により、入力された画像に対して画像復元処理の対象となる局所領域(画像の特定の範囲)を指定する。処理フレーム指定部33は、使用者の入力操作により、後述する処理フレーム数を指定する。
【0037】
記憶部40は、コンピュータの本体を構成する一要素であって、主記憶部41及び表示メモリ46を有する。主記憶部41は、入力インタフェース部20を介して取り込んだ画像データが記憶される取込画像メモリ42と、後述する画像復元処理を実行するための画像復元プログラムや各種パラメータが記憶された処理パラメータメモリ43とを有する。取込画像メモリ42は、少なくとも後述する画像復元処理において必要となる複数フレーム分の画像データを記憶可能な容量とする。
【0038】
表示メモリ46は、ディスプレイ等で構成された画像出力装置(図示省略)に表示するための画像データが記憶される領域であり、取り込んだ画像(劣化画像)を表示するために復元処理前の画像データ(取り込んだ画像データ)が記憶される取込画像メモリ47と、画像復元処理において一時的に必要な記憶領域である処理メモリ48と、復元された画像を表示するために復元処理後の画像データが記憶される復元画像メモリ49とを有する。
【0039】
また、図示は省略するが、記憶部40は、主記憶部41の容量不足を補う補助記憶部を設けてもよい。補助記憶部は、例えば、ハードディスク(HD)であってもよいし、CD−ROM、DVD、SSD(Solid State Drive)、フラッシュメモリなどの可搬性のあるものであってもよいし、また、それらの組み合わせであってもよい。
【0040】
また、本実施の形態では、処理パラメータメモリ43に予め画像復元プログラムが記憶されている場合について説明するが、画像復元プログラムは、記録メディアから入力インタフェース部20を介して記憶部40にインストールしたり、モデム及び入力インタフェース部20を介して外部から記憶部40にダウンロードしたりしてもよい。
【0041】
画像復元処理部50は、後述する画像復元処理を実行するための処理部であり、機能的には、
図2に示すように、位置変化量推定部51と、システム状態量算出部52と、観測量算出部53と、推定部54と、復元画像出力部55とを含んだ構成で表すことができる。
【0042】
位置変化量推定部51は、連続するフレームの画像間で対応する小領域の位置合わせを行い、小領域の位置変化量を推定する。
【0043】
システム状態量算出部52は、時刻nにおけるシステム状態量を、時刻(n−1)までに取得されたフレームから得られる時刻(n−1)のシステム状態量、位置変化量に基づく時刻(n−1)から時刻nへの状態遷移、及び時刻(n−1)のシステム状態量に含まれない時刻nのシステム状態量の要素を含む有色性駆動源を用いて算出する。
【0044】
観測量算出部53は、時刻nに取得されたフレームの画像情報、時刻nのシステム状態量、並びにぼけ及び雑音を示す情報を用いて、原画像情報にぼけ及び雑音が加わって劣化画像情報となる過程を観測量として算出する。
【0045】
推定部54は、算出されたシステム状態量及び観測量を有色性駆動源付カルマンアルゴリズムにより解いて、時刻nのシステム状態量の最適推定値を求める。
【0046】
復元画像出力部55は、推定部54により得られた時刻nのシステム状態量の最適推定値を、復元画像を示す画像情報として出力する。
【0047】
なお、画像復元処理部50の各部の処理については、後述する画像復元処理において詳述する。
【0048】
出力インタフェース部60は、画像復元処理部50で復元処理され、復元画像メモリ49に記憶された画像データ(復元画像)を、ディスプレイなどで構成された画像出力装置に出力可能なデータ形式に変換するなどの出力処理を行う。
【0049】
なお、復元画像の出力形態は画像出力装置に表示する場合に限定されず、他の画像出力装置に出力するようにしてもよい。例えば、画像出力装置として、プリンタ、記録メディア、モデムなどを用いることができる。記録メディア及びモデムは、画像入力装置としての記録メディア及びモデムとそれぞれ共用してもよい。なお、画像出力装置の種類は、画像復元装置10の用途に応じて適用可能なものを選択すればよい。
【0050】
また、図示は省略するが、出力インタフェース部60は、画像出力装置の種類に応じて別個独立に設けられる。上記のように、記録メディア及びモデムに関しては、通常、入力インタフェース部20と出力インタフェース部60とは一体化されている。出力インタフェース部60も、入力インタフェース部20と同様に、コンピュータ本体の内部に格納された内蔵カード(ボード)であってもよいし、内部インタフェース部を介して接続された外部設置型機器であってもよい。
【0051】
次に、第1の実施の形態における画像復元の理論となる劣化画像モデル及び状態空間モデルについて説明する。
【0052】
劣化した画像は、一般的に、原画像とぼけの点拡がり関数(PSF:Point Spread Function、以下単に「ぼけ関数」ともいう)の畳み込み(コンボリューション)に雑音を加えたモデルによって得られる。
【0053】
図3に、画像の劣化の過程を説明するための図を示す。例えば、
図3の左端に示す原画像にぼけが生じると、
図3の中央に示す画像となり、この画像にさらに雑音が加わると、
図3の右端に示す劣化画像となる。ぼけは、ある画素がその周囲の画素に影響することによって生じるが、雑音は、画素に関係なく発生する。カメラなどで撮影された画像の場合、ぼけは、手ぶれや焦点ずれなどに起因して生じ、雑音は、暗電流や熱雑音などに起因して不可避的に発生する。従って、第1の実施の形態における劣化画像の復元は、劣化画像からぼけ及び雑音を取り除くことである。
【0054】
第1の実施の形態では、AR係数を用いることなく、クリアな画像情報(原画像情報)のみを用いてシステム状態量を算出し、かつ、劣化画像情報、クリアな画像情報、ぼけ情報、及び雑音を用いて観測量を算出した状態空間モデルにおいて、有色性駆動源付カルマンアルゴリズムにより、システム状態量の最適推定値を求めることにより、劣化画像を復元する。
【0055】
上記の状態空間モデルは、下記(1)式の状態方程式と、下記(2)式の観測方程式とで表される。ただし、状態方程式において、ベクトルx
p1(n)は時刻nの状態ベクトル(原画像のシステム状態量)、行列Φ
p1は状態遷移行列、ベクトルδ
p1(n)は駆動源ベクトルである。また、観測方程式において、ベクトルy
p1(n)は観測ベクトル(観測画像情報、劣化画像情報)、行列M
p1は観測行列、ベクトルε
p1(n)は観測雑音ベクトルである。また、添え字p1は第1の実施の形態に係るものであることを表す。
【0057】
(1)式の状態方程式は、観測対象のシステムを状態空間モデルで記述したものであり、内部状態つまり状態変数(ここでは、状態ベクトルx
p1(n−1)、x
p1(n))の時間に対する生成過程を表している。また、(2)式の観測方程式は、何らかの観測装置を通じて観測する過程を記述したものであり、観測結果(ここでは、観測ベクトルy
p1(n))が、被観測量つまり入力(ここでは、状態ベクトルx
p1(n))に依存して生成される様子を示している。
【0058】
次に、
図4を参照して、第1の実施の形態に係る画像復元装置10の作用について説明する。なお、
図4に示す画像復元処理ルーチンは、記憶部40に制御プログラムとして格納されており、図示しないCPUによって実行される。
【0059】
ステップ100で、画像入力装置から入力インタフェース部20を介して、所定フレームの画像情報を取り込み、主記憶部41の取込画像メモリ42に記憶する。また、表示メモリ46の取込画像メモリ47にも記憶し、出力インタフェース部60を介して、画像出力装置に取込画像を表示する。なお、以下では、時刻nのフレームの画像情報をLフレームと表記する。また、ここでは、時刻nの開始値を「0」とし、n≧(k−1)となる時刻のフレームから、以下の処理を行う。kは後述する処理フレーム数であり、本実施の形態では「3」である。ここでは、n=2から以下の処理を行う。
【0060】
次に、ステップ102で、位置変化量推定部51が、取込画像メモリ42に記憶されている連続するフレームの画像情報を2つずつ読み出し、各フレームの画像間の位置合わせを行い、連続するフレーム間における画像内の物体の位置変化量を推定する。画像の位置合わせには、ブロックマッチングやLucas−Kanade法といった既存の手法を用いることができる。本実施の形態では、
図5に示すようなブロックマッチングを用いる場合について説明する。
【0061】
ブロックマッチングとは、時刻(n−1)の(L−1)フレームで設定した探索範囲内において、時刻nのLフレームで類似する領域を探す手法であり、時刻(n−1)の(L−1)フレームの小領域の画素値と時刻nのLフレームで取り出した小領域の画素値との誤差の総和が最小となる領域を求め、対応する位置を探す手法である。画素値の誤差の総和を取る手段として、下記(3)式に示す各画素値の2乗距離の総和を取るSSD(Sum of Squared Distance)や、下記(4)式に示す各画素値の絶対値距離の総和を取るSAD(Sum of Absolute Distance)が挙げられる(非特許文献2「藤原 慎矢、田口 亮、“ブロックマッチング動き推定に基づくフレームレート向上に関する一手法”、IEICE Technical Report SIS2006-15, Jun. 2006」参照)。なお、(3)式及び(4)式におけるx
i,jは画素(i,j)の画素値である。
【0063】
例えば、画像内の物体が左から右へ水平方向に移動していた場合には、時刻nのLフレームにおいて、時刻(n−1)の(L−1)フレームに設定した小領域の位置より右方向にずれた位置から取り出した小領域が、時刻(n−1)の(L−1)フレームの小領域とマッチングすることになる。このときのずれ量を時刻(n−1)の(L−1)フレームと時刻nのLフレーム間の位置変化量として推定する。以下、画像内の物体が左から右へ水平方向に移動していた場合を例に説明する。
【0064】
次に、ステップ104で、位置変化量推定部51が、画像復元処理の処理単位となる局所領域のサイズを設定する。局所領域のサイズは、上記ステップ102で推定された位置変化量に応じて決定される。局所領域のサイズの設定は、ユーザにより操作部30の領域指定部32から指定された値を設定するようにしてもよいし、位置変化量と局所領域のサイズとの関係を予めテーブル等に定めておいて、推定された位置変化量に基づいて自動的に設定するようにしてもよい。
【0065】
また、位置変化量推定部51が、後述する状態方程式における時刻(n−1)のシステム状態量として用いるフレーム数を設定する。フレーム数の設定は、ユーザにより操作部30の処理フレーム指定部33から指定された値を設定するようにしてもよいし、予め定めた値を設定するようにしてもよい。
【0066】
ここでは、ステップ104において、局所領域のサイズp×q=3×3画素、フレーム数3と設定されたものとする。この場合、
図6に示すように、左から右へ流れている、時刻(n−2)の(L−2)フレーム、時刻(n−1)の(L−1)フレーム、時刻nのLフレームという時系列の3フレーム分の画像情報が処理対象となる。
【0067】
次に、ステップ106で、位置変化量推定部51が、上記ステップ104で設定されたサイズの局所領域を各フレームに対して同一の領域に設定する。
【0068】
次に、ステップ108で、システム状態量算出部52が、システム状態量算出処理を実行する。ここで、
図7を参照して、システム状態量算出処理ルーチンについて説明する。
【0069】
ステップ1080で、システム状態量算出部52が、時刻(n−2)の(L−2)フレーム及び時刻(n−1)の(L−1)フレームに設定された局所領域内の各画素のシステム状態量を要素とする時刻(n−1)の状態ベクトルと、時刻nのLフレームに設定された局所領域内の各画素のシステム状態量を要素とする時刻nの状態ベクトルとを導出する。以下、各時刻の状態ベクトルの導出方法を具体的に説明する。
【0070】
まず、時刻(n−1)の状態ベクトルは、
図8に示すように、時刻(n−1)の(L―1)フレーム及び時刻(n−2)の(L―2)フレームにおけるp×qサイズの局所領域内の各画素のシステム状態量を組み合わせて構成される。すなわち、時刻(n−1)のシステム状態ベクトルは18次元のベクトルx
p1(n−1)=[x
1,1L−1(n−1) x
2,1L−1(n−1) ・・・ x
3,3L−1(n−1) x
1,1L−2(n−1) x
2,1L−2(n−1) ・・・ x
3,3L−2(n−1)]
Tとなる。次に、時刻nの状態ベクトルは、
図9に示すように、時刻nのLフレームにおけるp×qサイズの局所領域内の各画素のシステム状態量を組み合わせて構成される。すなわち、時刻nの状態ベクトルは9次元のベクトルx
p1(n)=[x
1,1L(n) x
2,1L(n) ・・・ x
3,3L(n)]
Tとなる。ただし、Tはベクトルの転置を表す。また、状態ベクトルx
p1(n−1)及びx
p1(n)のサイズ(次元数)は設定されたフレーム数及び局所領域p×qのサイズによって変化する。
【0071】
次に、ステップ1082で、システム状態量算出部52が、上記ステップ102で推定した位置変化量、及び時刻nのシステム状態量と時刻(n−1)のシステム状態量との相関に基づく忘却係数を用いて状態遷移行列Φ
p1を構成する。
【0072】
具体的には、まず、画像内の物体が、時刻(n−2)から時刻(n−1)間、及び時刻(n−1)から時刻n間に、水平に左から右へ1画素ずつ移動していたとする。この場合、各フレームの局所領域において、
図10に示すように、同一のシステム状態量を持つ領域が存在する。ここで、同一のシステム状態量を持つ領域は、上記ステップ102で推定された位置変化量から判断できる。同一のシステム状態量を持つ領域は、
図10において網掛けされた画素領域である。
図11に、
図10の各フレームの局所領域部分を抜き出した図を示す。
図11に示すように、時刻(n−2)の(L−2)フレームのシステム状態量x
1,1L−2,x
2,1L−2,x
3,1L−2と、時刻(n−1)の(L−1)フレームのシステム状態量x
1,2L−1,x
2,2L−1,x
3,2L−1と、時刻nのLフレームのシステム状態量x
1,3L,x
2,3L,x
3,3Lとが同一のシステム状態量を持つ領域となる。また、時刻(n−1)の(L−1)フレームのシステム状態量x
1,1L−1,x
2,1L−1,x
3,1L−1と、時刻nのLフレームのシステム状態量x
1,2L,x
2,2L,x
3,2Lとも同一のシステム状態量を持つ領域となる。
【0073】
また、時刻(n−1)のシステム状態量と時刻nのシステム状態量とは相関を持ち、
図12に示すように、時刻(n−1)の状態ベクトルを構成する時刻(n−2)の(L−2)フレーム及び時刻(n−1)の(L−1)フレームにおいて、時刻nのLフレームとの時刻が近いほど強い相関を示し、遠いほど弱い相関を示す。そこで、時刻nのLフレームからフレームが遠ざかるにつれて相関が低くなるようなパラメータ(忘却係数ρ
i,j)を設定する。一般的に忘却係数は時刻nのLフレームに近いほど相関が強まるように1に近い値となり、時刻nのLフレームから遠ざかるほど相関が弱まるように0に近い値となるように構成する。ここで、忘却係数は推定して導出してもよいし、既知として扱ってもよい。
【0074】
そして、上述した同一のシステム状態量を持つ領域が時刻(n−1)から時刻nへの変化によって一致し、かつ忘却係数を導入した状態遷移行列Φ
p1を構成する。
【0075】
次に、ステップ1084で、システム状態量算出部52が、駆動源ベクトルを構成する。
図11の時刻nのLフレームのシステム状態量x
1,1L,x
2,1L,x
3,1Lは、時刻(n−1)の情報からは構成することができないため、本実施の形態では、駆動源と呼ばれる領域にシステム状態量x
1,1L,x
2,1L,x
3,1Lを用いることで、駆動源ベクトルδ
p1(n)を構成する。これは、物体がフレーム内を位置変化した際に、時刻(n−1)のシステム状態量((時刻(n−2)のL−2)フレーム及び時刻(n−1)の(L−1)フレーム)では表せられない領域のシステム状態量を、時刻nのシステム状態量で構成することを表しており、原画像の情報が含まれることから、駆動源ベクトルδ
p1(n)は有色性駆動源となる。
【0076】
次に、ステップ1086で、システム状態量算出部52が、上記ステップ1080で導出した時刻(n−1)及び時刻nの状態ベクトルx
p1(n−1)及びx
p1(n)、上記ステップ1082で構成した状態遷移行列Φ
p1、及び上記ステップ1084で構成した駆動源ベクトルδ
p1(n)に基づいて、システム状態量を算出する。算出されたシステム状態量(状態方程式)の一例と、各ステップで構成された各情報との関係を
図13に示す。このように、本実施の形態における状態方程式は、クリアな画像からなる状態ベクトルx
p1(n)と、0及び忘却係数ρ
i,jからなる状態遷移行列Φ
p1と、有色信号であるクリアな画像からなる駆動源ベクトルδ
p1(n)とから構成されるため、所望のシステム状態量(原画像の画素情報)、つまり、クリアな画像情報(原画像情報)のみで構成されることになる。
【0077】
システム状態量算出処理ルーチンが終了すると、画像復元処理ルーチンにリターンして、ステップ110へ移行し、観測量算出部53が、
図3に示した劣化画像モデルに基づいて、
図14に示すように、原画像情報にぼけ及び雑音が加わって劣化画像情報となる過程を、観測量(観測方程式)として算出する。観測ベクトルy
p1(n)は、劣化画像情報を含む観測量として時刻nのLフレームの局所領域内の各画素の画素情報を要素とする9次元のベクトルである。また、観測行列M
p1は、画像の劣化モデルにおけるぼけの点拡がり関数(PSF)に対応したh
i,jを要素とする行列である。各要素h
i,j(i,jはhの座標であり、i,j=−1,0,1)は、既知であって、予めデータ化して適切に設定されている。また、観測雑音ベクトルε
p1(n)は、局所領域に含まれる画素に対応する観測雑音を要素とする9次元のベクトルである。
【0078】
次に、ステップ112で、推定部54が、上記ステップ108で算出されたシステム状態量(状態方程式)及び上記ステップ110で算出された観測量(観測方程式)により、下記に示す有色性駆動源付カルマンアルゴリズムを導出する。
【0080】
上記のアルゴリズムは、初期設定の過程[Initialization]と反復の過程[Iteration]とに大別され、反復の過程では、1〜5の手順を逐次繰り返す。以下、
図15を参照して、有色性駆動源付カルマンアルゴリズムの詳細について説明する。
【0081】
ステップ1120で、カルマンアルゴリズムの処理が初回か否かを判定する。初回の場合には、ステップ1122へ移行し、2回目以降の場合には、ステップ1124へ移行する。
【0082】
ステップ1122では、初期設定を行う。具体的には、状態ベクトル、つまり、システム状態量である所望信号(原画像信号)ベクトルの最適推定値(以下「所望信号の最適推定値ベクトル」という)の初期値x^
p1(0|0)、所望信号ベクトルの推定誤差(以下「所望信号の推定誤差ベクトル」という)の相関行列の初期値P
p1(0|0)、駆動源ベクトルの共分散の値R
δp1(n)[i,j]、及び観測雑音ベクトルの共分散の値R
εp1(n)[i,j]を、上述の初期設定の過程[Initialization]に示した初期状態に設定する。
【0083】
ただし、ベクトル0
2K^2は、2K
2次元の零ベクトルであり、行列I
2K^2は、2K
2次の単位行列である。また、Kは、K×Kの大きさの局所領域のその次元である。また、E[δ
p1(n),δ
p1T(n)]は、駆動源ベクトルδ
p1(n)の自己相関行列である。E[ε
p1(n)ε
p1T(n)][i,j]は、観測雑音ベクトルε
p1の自己相関行列であり、ここではσ
2ε[i,j]と等しく、既知と仮定している。ここでいう「既知」とは、別の任意の方法(アルゴリズム)で求められて与えられることを意味する。もし観測雑音ベクトルε
p1(n)が白色雑音であり零平均であれば、σ
2εは所定のサンプル数の分散で与えられる。なお、ベクトル及び行列の要素を示す場合、ベクトルa(n)のi番目の要素をa(n)[i]と表記し、行列A(n)のi行j列の要素をA(n)[i,j]と表記することとする。
【0084】
次に、ステップ1124で、上記ステップ1082で定義した状態空間モデルにおける状態遷移行列Φ
p1、上記ステップ1122で設定した所望信号の推定誤差ベクトルの相関行列の初期値P
p1(0|0)(n=2の場合)、または1時刻前に後述するステップ1130で更新された相関行列P
p1(n|n)(n>2の場合)、及び駆動源ベクトルの共分散の値R
δp1(n)[i,j]を用いて、時刻(n−1)までの情報により時刻nの状態ベクトルを推定した場合の誤差である相関行列P
p1(n|n−1)を計算する(上述の反復の過程[Iteration]の手順1)。
【0085】
次に、ステップ1126で、上記ステップ1124で計算した相関行列P
p1(n|n−1)、上記ステップ110で定義した状態空間モデルにおける観測行列M
p1、及び観測雑音ベクトルの共分散の値R
εp1(n)[i,j]を用いて、カルマンゲイン行列K
p1(n)を計算する(同手順2)。
【0086】
次に、ステップ1128で、状態遷移行列Φ
p1、及び上記ステップ1122で設定した所望信号の最適推定値ベクトルの初期値x^
p1(0|0)(n=2の場合)、または1時刻前に本ステップで得られた所望信号の最適推定値ベクトルx^
p1(n|n)(n>2の場合)を用いて、時刻(n−1)までの情報による時刻nでの所望信号の最適推定値ベクトルx^
p1(n|n−1)を計算する(同手順3)。そして、計算した所望信号の最適推定値ベクトルx^
p1(n|n−1)、上記ステップ1126で計算したカルマンゲイン行列K
p1(n)、観測ベクトルy
p1(n)、及び観測行列M
p1を用いて、時刻nまでの情報によるその時刻での所望信号の最適推定値ベクトルx^
p1(n|n)を計算する(同手順4)。
【0087】
次に、ステップ1130で、単位行列I、カルマンゲイン行列K
p1(n)、観測行列M
p1、及び上記ステップ1124で計算された相関行列P
p1(n|n−1)を用いて、時刻nまでの情報によるその時刻での相関行列P
p2(n|n)を更新する。次に、ステップ1132で、上記ステップ1128で計算された所望信号の最適推定値ベクトルx^
p1(n|n)を、現在設定されている局所領域の処理結果として一旦復元画像メモリ49に記憶する。
【0088】
有色性駆動源付カルマンアルゴリズムが終了すると、画像復元処理へリターンして、ステップ114へ移行し、現時刻のフレームの全領域について画像復元処理が終了したか否かを判定する。全領域について処理が終了していない場合には、ステップ106へ戻って、次の局所領域を設定して処理を繰り返す。全領域について処理が終了した場合には、ステップ116へ移行して、上記ステップ100で取り込んだ全フレームについて処理が終了したか否かを判定する。全フレームについて処理が終了していない場合には、ステップ118へ移行してnを1インクリメントして、ステップ102へ戻って、次のフレームについて処理を繰り返す。全フレームについて処理が終了した場合には、ステップ120へ移行する。
【0089】
ステップ120では、復元画像出力部55が、復元画像メモリ49に記憶された画像データを、出力インタフェース部60を介して出力して、画像出力装置に復元画像を表示し、画像復元処理を終了する。
【0090】
以上説明したように、第1の実施の形態の画像復元装置10によれば、時系列の複数フレームの画像を用いて、画像内の物体の位置変化量を推定し、この位置変化量に基づく状態遷移行列及び駆動源ベクトルを導入したシステム状態量を算出することにより、AR係数を推定することなく、物体の急激な変化にも対応して、劣化画像を精度良く復元することができる。
【0091】
また、状態遷移行列を構成する際に、時刻nからの時間間隔に応じた忘却係数を導入することで、状態遷移行列をより適切に表現することができ、より精度良く画像を復元することができる。
【0092】
また、局所領域毎に処理していくことで、処理する行列のサイズを小さくして演算量を抑えることができ、並列処理を適用することで、画像全体の復元処理を高速に実行することができる。
【0093】
なお、第1の実施の形態では、画像内の物体が水平方向に移動する場合を例に説明したがこれに限定されない。位置変化量推定部51で推定された位置の変化に応じたシステム状態量を算出すればよい。例えば、
図16に示すように、画像内の物体が上方へ移動する場合でも、システム状態量が重なる領域を考慮した状態遷移行列Φ
p1及び駆動源ベクトルδ
p1を構成するとよい。
【0094】
また、第1の実施の形態では、忘却係数として時刻変化における相関値の変化に基づく値を用いる場合について説明したが、時刻変化によらないフレーム間の相関値に基づく値を用いるようにしてもよい。
【0095】
次に、第1の実施の形態に係る画像復元装置10の効果を実証するためのシミュレーションについて説明する。本シミュレーションでは、第1の実施の形態に係る画像復元装置10の画像復元能力を評価するために、(1)視覚評価、(2)客観評価、及び(3)主観評価を行った。視覚評価は、原画像と復元画像とを視覚により比較した評価である。客観評価は、数値による評価である。主観評価は、聞き取り調査である。以下、これらを順に説明する。
【0096】
まず、シミュレーションの条件について説明する。本シミュレーションでは、原画像として、
図17に示す2つの画像(「Girl」及び「Building」)を用いた。いずれも階調8ビットのグレースケール画像で、解像度は256×256画素である。また、原画像に加えるぼけに対応する点拡がり関数(PSF)のモデルとして、2次元ガウス関数を用い、雑音として、加法性白色ガウス雑音を用いた。また、信号対雑音比SNRを20dBとした。また、画像復元に用いる複数フレームの画像として、水平方向に3画素移動させた画像を用いた。また、同一の上記シミュレーション条件の下で、第1の実施の形態における手法(発明手法1)と非特許文献1による手法(従来手法1)とを比較した。
【0097】
(1)視覚評価
画像「Girl」についてのシミュレーション結果を
図18及び
図19に、画像「Building」についてのシミュレーション結果を
図20及び
図21に示す。なお、
図19は
図18の一部を拡大したもの、
図21は
図20の一部を拡大したものである。
【0098】
図19に示すように、従来手法1は、輪郭情報は復元できているものの、雑音による劣化が強調されている。一方、発明手法1は、輪郭情報を維持しつつ鮮明な画像を復元することができている。また、
図21に示すように、従来手法1は、樹木や建物の線が劣化している。一方、発明手法1は、ぼけが除去され、鮮明な画像を復元できている。このように、発明手法1の有効性を確認することができる。
【0099】
(2)客観評価(数値による評価)
下表1に、原画像に対するシミュレーション結果(客観評価)を示す。
【0101】
ここでは、従来手法1及び発明手法1の画像復元能力を数値により評価するため、下記(5)式に示すPSNR[dB]を用いて画像復元能力を評価した。なお、PSNR(Peak Signal-to-Noise Ratio)は、画像の信号に対する雑音の比であり、数値が大きいほど劣化(ぼけや雑音など)が少なく画像として良好であるといえる。
【0103】
表1に示すように、「Girl」及び「Building」のどちらの画像についても、発明手法1は、従来手法1よりもPSNRの数値が大きいことが確認できる。これにより、発明手法1は、客観評価の点でも、従来手法1に比べて画像復元能力が高いことがわかる。
【0104】
(3)主観評価(聞き取り調査)
図22に、原画像に対するシミュレーション結果(主観評価)を示す。
【0105】
ここでは、従来手法1及び発明手法1の画像復元能力を評価するために、聞き取り調査による主観評価を行った。画像復元の性能評価は、下表2に示すように、MOS(Mean Opinion Score)評価基準に基づく5段階の評価値を用いた聞き取り調査により行った。50人(男性25人、女性25人)の聴取者が画像復元により得られた画像(
図18〜
図21参照)を評価した。各々の聴取者は、評価値「1」から「5」を決定する。評価値「5」が最良である。
【0107】
図22に示すように、「Girl」及び「Building」のどちらの画像についても、発明手法1は、従来手法1よりも高いMOS評価値を得ていることが確認できる。これにより、発明手法1は、主観評価の点でも、従来手法1に比べて画像復元能力が高いことがわかる。
【0108】
以上の実験結果により、第1の実施の形態に係る画像復元装置10における手法(発明手法1)は、従来手法1に比べて高い画像復元能力を発揮していることがわかる。
【0109】
<第2の実施の形態>
次に、第2の実施の形態について説明する。第2の実施の形態では、両面に画像がプリントされた原稿をスキャンした際の裏写りを除去する画像復元装置に本発明を適用した場合について説明する。
【0110】
第2の実施の形態に係る画像復元装置210の構成について、第1の実施の形態に係る画像復元装置10と異なる部分について説明し、他の部分については第1の実施の形態と同一符号を付して詳細な説明を省略する。
【0111】
図23に示すように、第2の実施の形態における画像復元処理部250は、位置変化量設定部251と、システム状態量算出部252と、観測量算出部253と、推定部54と、復元画像出力部55とを含んだ構成で表すことができる。
【0112】
位置変化量設定部251は、裏面画像のフレーム及び表面画像のフレームを時刻(n−1)までに取得されたフレームとし、表面画像の原画像を時刻nに取得されたフレームとし、時刻(n−1)から時刻n間の表面画像内の物体の位置変化量を想定した局所領域を設定する。
【0113】
システム状態量算出部252は、時刻nにおけるシステム状態量を、時刻(n−1)までに取得されたフレームから得られる時刻(n−1)のシステム状態量、位置変化量に基づく時刻(n−1)から時刻nへの状態遷移、及び時刻(n−1)のシステム状態量に含まれない時刻nのシステム状態量の要素を含む有色性駆動源を用いて構成する。
【0114】
観測量算出部253は、表面画像のフレームの画像情報、時刻nのシステム状態量、及び裏写りに対応する劣化情報を用いて、表面画像の原画像情報に裏面画像の裏写りが生じて劣化画像情報となる過程を観測量として算出する。
【0115】
次に、第2の実施の形態における画像復元の理論となる劣化画像モデル及び状態空間モデルについて説明する。
【0116】
図24に、画像の劣化の過程を説明するための図を示す。例えば、両面に画像がプリントされた原稿をスキャンした場合、表面画像の原画像情報に裏面画像の画像情報が裏写りして、
図24の右端に示す劣化画像となる。裏写りは、原画像(表面画像)をスキャンした画像に、裏面画像が透過し、かつぼけが生じた状態の画像が重畳された画像である。なお、裏面画像が透過する透過率は、原稿の厚さ、画像の濃度、スキャン時の照明強度等に依存する。従って、第2の実施の形態にける劣化画像の復元は、劣化画像から裏写りした裏面画像情報を取り除くことである。
【0117】
第2の実施の形態においても、AR係数を用いることなくシステム状態量を算出し、かつ裏写りによる画像の劣化過程に基づいて観測量を算出した状態空間モデルにおいて、有色性駆動源付カルマンアルゴリズムにより、システム状態量の最適推定値を求めることにより、劣化画像を復元する。
【0118】
上記の状態空間モデルは、下記(6)式の状態方程式と、下記(7)式の観測方程式とで表される。ただし、状態方程式において、ベクトルx
p2f(n)は時刻nの状態ベクトル(表面画像の原画像のシステム状態量)、行列Φ
p2は状態遷移行列、x
p2(n−1)は時刻(n−1)の状態ベクトル(表面画像及び裏面画像のシステム状態量)、ベクトルδ
p2(n)は駆動源ベクトルである。また、観測方程式において、ベクトルy
p2f(n)は観測ベクトル(観測画像情報、劣化画像情報)、行列M
p2は観測行列、HはぼけのPSF、ベクトルx
p2b(n)は時刻nの状態ベクトル(裏面画像のシステム状態量)である。また、添え字p2は第2の実施の形態に係るものであることを表す。
【0120】
次に、
図25を参照して、第2の実施の形態に係る画像復元装置210の作用について説明する。なお、
図25に示す画像復元処理ルーチンは、記憶部40に制御プログラムとして格納されており、図示しないCPUによって実行される。また、第1の実施の形態における画像復元処理ルーチンと同様の処理については、同一符号を付して詳細な説明を省略する。
【0121】
ステップ200で、画像入力装置であるスキャナにより、両面に画像がプリントされた原稿から読み取られた表面画像及び裏面画像の画像情報を、入力インタフェース部20を介して取り込み、主記憶部41の取込画像メモリ42に記憶する。以下では、第1の実施の形態に対応させて、表面画像の原画像をLフレーム、読み取った表面画像を(L−1)フレーム、及び裏面画像を(L−2)フレームとする。なお、(L−2)フレームは、読み取られた裏面画像を左右(または上下)反転させた画像とする。
【0122】
次に、ステップ202で、位置変化量設定部251が、
図26に示すように、予め定められたサイズ(例えば、p×q=3×3画素)の局所領域を、裏面画像の(L−2)フレーム及び表面画像の(L−1)フレームの同一の領域に設定する。また、表面画像の原画像のLフレームにおいて、裏面画像の(L−2)フレーム及び表面画像の(L−1)フレームに局所領域を設定した位置から予め定めた位置変化量分ずらした領域に局所領域を設定する。この位置変化量は、第1の実施の形態と同様に、時間変化における表面画像内の物体の位置変化を示すものであり、ここでは、位置変化量を右方向へ1画素分とする。
【0123】
次に、ステップ204で、システム状態量算出部252が、システム状態量算出処理を実行する。システム状態量算出処理は、第1の実施の形態におけるシステム状態量算出処理と処理の流れ自体は同様であり、各処理の詳細が異なるため、第1の実施の形態で用いた
図7を参照して、システム状態量算出処理ルーチンについて説明する。
【0124】
ステップ1080で、システム状態量算出部252が、裏面画像の(L−2)フレーム及び表面画像の(L−1)フレームに設定された局所領域内の各画素のシステム状態量を要素とする時刻(n−1)の状態ベクトルと、表面画像の原画像のLフレームに設定された局所領域内の各画素のシステム状態量を要素とする時刻nの状態ベクトルとを導出する。以下、各時刻の状態ベクトルの導出方法を具体的に説明する。
【0125】
まず、時刻(n−1)の状態ベクトルは、
図27に示すように、表面画像の(L―1)フレーム及び裏面画像の(L―2)フレームにおけるp×qサイズの局所領域内の各画素のシステム状態量を組み合わせて構成される。すなわち、時刻(n−1)の状態ベクトルは18次元のベクトルx
p2(n−1)=[x
1,1L−1(n−1) x
2,1L−1(n−1) ・・・ x
3,3L−1(n−1) x
1,1L−2(n−1) x
2,1L−2(n−1) ・・・ x
3,3L−2(n−1)]
Tとなる。次に、時刻nの状態ベクトルは、
図28に示すように、表面画像の原画像のLフレームにおけるp×qサイズの局所領域内の各画素のシステム状態量を組み合わせて構成される。すなわち、時刻nの状態ベクトルは9次元のベクトルx
p2(n)=[x
1,1L(n) x
2,1L(n) ・・・ x
3,3L(n)]
Tとなる。ただし、Tはベクトルの転置を表す。また、状態ベクトルx
p2(n−1)及びx
p2(n)のサイズ(次元数)は設定された局所領域p×qのサイズによって変化する。
【0126】
次に、ステップ1082で、システム状態量算出部252が、位置変化量(局所領域の設定位置のずれ量)、及び時刻nのシステム状態量と時刻(n−1)のシステム状態量との相関に基づく忘却係数を用いて状態遷移行列Φ
p2を構成する。
【0127】
図29に示すように、時刻(n−1)の状態ベクトルを構成する裏面画像の(L−2)フレーム及び表面画像の(L−1)フレームにおいて、表面画像の(L−1)フレームの方が表面画像の原画像であるLフレームと強い相関を示し、裏面画像の(L−2)フレームの方が弱い相関を示す。そこで、表面画像の(L−1)フレームの要素に対しては相関が高く、裏面画像の(L−2)フレームの要素に対しては相関が低くなるようなパラメータ(忘却係数ρ
i,j)を設定する。ここで、忘却係数は推定して導出してもよいし、既知として扱ってもよい。
【0128】
そして、位置変化量に基づいて得られる同一のシステム状態量を持つ領域が時刻(n−1)から時刻nへの変化によって一致し、かつ忘却係数を導入した状態遷移行列Φ
p2を構成する。
【0129】
次に、ステップ1084で、システム状態量算出部252が、駆動源ベクトルを構成する。ここでは、表面画像の(L−1)フレーム及び裏面画像の(L−2)フレームに対して、表面画像の原画像のLフレームの局所領域を右方向にずらしているため、表面画像の原画像のLフレームのシステム状態量x
1,4L,x
2,4L,x
3,4Lは、時刻(n−1)の情報からは構成することができないため、本実施の形態では、駆動源と呼ばれる領域にシステム状態量x
1,4L,x
2,4L,x
3,4Lを用いることで、駆動源ベクトルδ
p2(n)を構成する。
【0130】
次に、ステップ1086で、システム状態量算出部252が、上記ステップ1080で構成した時刻(n−1)及び時刻nの状態ベクトルx
p2(n−1)及びx
p2(n)、上記ステップ1082で構成した状態遷移行列Φ
p2、及び上記ステップ1084で構成した駆動源ベクトルδ
p2(n)に基づいて、システム状態量(状態方程式)を算出する。
【0131】
システム状態量算出処理ルーチンが終了すると、画像復元処理ルーチンにリターンして、ステップ206へ移行し、観測量算出部253が、
図24に示した劣化画像モデルに基づいて、
図30に示すように、表面画像の原画像に裏面画像の裏写りが加わって劣化画像となる過程を、観測量(観測方程式)として算出する。観測ベクトルy
p2f(n)は、劣化画像情報を含む観測量として表面画像のフレームの局所領域内の各画素の画素情報を要素とする9次元のベクトルである。また、観測行列M
p2は、観測方程式を満たすための単位行列である。ベクトルx
p2f(n)は、表面画像の原画像における局所領域内の各画素のシステム状態量を要素とするベクトルである。sは裏写りに関する透過率であり、ユーザにより操作部30のパラメータ設定部31から設定された値を用いる。透過率は、上述のように、原稿の厚さ、画像の濃度、スキャン時の照明の強度等に応じた値を設定しておく。行列Hは、画像の劣化モデルにおけるぼけの点拡がり関数(PSF)に対応したh
i,jを要素とする行列である。ベクトルx
p2b(n)は、裏面画像の局所領域内の各画素のシステム状態量を要素とする9次元のベクトルである。
【0132】
次に、ステップ112で、推定部54が、上記ステップ204で算出されたシステム状態量(状態方程式)及び上記ステップ206で算出された観測量(観測方程式)により、第1の実施の形態と同様に、有色性駆動源付カルマンアルゴリズムを実行する。以下、ステップ114〜120を実行して、画像復元処理を終了する。
【0133】
以上説明したように、第2の実施の形態の画像復元装置210によれば、表面画像及び裏面画像という複数フレームの画像を用いて、時刻変化における物体の位置変化を考慮した状態遷移行列及び駆動源ベクトルを導入したシステム状態量を算出することにより、AR係数を推定することなく、裏写りを精度良く除去することができる。
【0134】
また、状態遷移行列を構成する際に、表面画像の原画像と表面画像及び裏面画像との相関を示す忘却係数を導入することで、状態遷移行列をより適切に表現することができ、より精度良く画像を復元することができる。
【0135】
また、局所領域毎に処理していくことで、処理する行列のサイズを小さくして演算量を抑えることができ、並列処理を適用することで、画像全体の復元処理を高速に実行することができる。
【0136】
次に、第2の実施の形態に係る画像復元装置210の効果を実証するためのシミュレーションについて説明する。本シミュレーションでも、発明手法1に対するシミュレーションと同様に、(1)視覚評価、(2)客観評価、及び(3)主観評価を行った。
【0137】
まず、シミュレーションの条件について説明する。本シミュレーションでは、
図31に示すように、表面の原画像として2つの画像(「Image1」及び「Lenne」)を用い、裏面画像として1つの画像(「Image2」)を用いた。「Image1」及び「Image2」は階調8ビットのグレースケール画像、「Lenna」は階調8ビットのカラースケール画像である。また、いずれの画像も解像度は256×256画素である。その他の条件は発明手法1に対するシミュレーションと同じである。また、同一の上記シミュレーション条件の下で、第2の実施の形態における手法(発明手法2)と非特許文献3(杉山 賢二、「実践 映像信号処理」、コロナ社、2008.)による手法(従来手法2)、及び特許文献4(特開2009−224873号公報)による手法(従来手法3)とを比較した。
【0138】
なお、従来手法2は、ガンマ補正を用いた補正手法であるが、表面画像の薄い色に対しても補正処理を行ってしまうため、カラー画像などに対しては精度良く裏写りを除去できない場合がある、という問題がある。
【0139】
また、従来手法3は、YUV(Y:輝度、U,V:色素)情報を用いてカラー画像を輝度変化してグレースケール化してエッジ検出して、所定の閾値により2値化処理し、2値化情報に基づいて色情報を埋め込む手法である。文字等のエッジが高い部分は精度良く裏写りを除去することができるが、広範囲に渡る裏写りに対しては処理が困難である、という問題がある。また、2値化処理の際の閾値を最適な一意に決定することは困難である、という問題もある。
【0140】
(1)視覚評価
画像「Image1」についてのシミュレーション結果を
図32及び
図33に、画像「Lunna」についてのシミュレーション結果を
図34及び
図35に示す。なお、
図33は
図32の一部を拡大したもの、
図35は
図34の一部を拡大したものである。
【0141】
図33に示すように、従来手法2は、裏写りが完全に除去できておらず不完全であり、従来手法3は、裏写りは除去できているものの、原画像の情報がわずかに劣化している。一方、発明手法2は、原画像の情報を維持しつつ、裏写りを精度良く除去できている。また、
図35に示すように、従来手法2は、画像全体が白くなり、また裏写りの除去も不完全であり、従来手法3も、裏写りが完全には除去できておらず不鮮明である。一方、発明手法2は、裏写りが除去できており鮮明な画像となっている。このように、発明手法2の有効性を確認することができる。
【0142】
(2)客観評価(数値による評価)
下表3に、(5)式に示したPSNRを用いたシミュレーション結果(客観評価)を示す。
【0144】
表3に示すように、「Image1」及び「Lenna」のどちらの画像についても、発明手法2は、従来手法2及び従来手法3よりもPSNRの数値が大きいことが確認できる。これにより、発明手法2は、客観評価の点でも、従来手法2及び従来手法3に比べて画像復元能力が高いことがわかる。
【0145】
(3)主観評価(聞き取り調査)
図36に、表2に示したMOS評価基準を用いたシミュレーション結果(主観評価)を示す。
【0146】
図36に示すように、「Image1」及び「Lenna」のどちらの画像についても、発明手法2は、従来手法2及び従来手法3よりも高いMOS評価値を得ていることが確認できる。これにより、発明手法2は、主観評価の点でも、従来手法2及び従来手法3に比べて画像復元能力が高いことがわかる。
【0147】
以上の実験結果により、第2の実施の形態に係る画像復元装置10における手法(発明手法2)は、従来手法2及び従来手法3に比べて高い画像復元能力を発揮していることがわかる。
【0148】
なお、上記実施の形態では、CPUによりプログラムを実行することで各部の処理を実行する場合について説明したが、各部をハードウエアにより構成するようにしてもよい。