(57)【特許請求の範囲】
【請求項１】
原画像情報に不必要な情報が混在した劣化画像情報のみから前記原画像情報を推定する画像復元装置であって、
時刻ｎのみの劣化画像情報に対して、時刻ｎまたは時刻ｎ＋１までの情報により前記原画像情報を含む時刻ｎ＋１でのシステムの状態量を推定した場合の推定誤差の相関値を算出する相関演算部と、
時刻ｎのみの劣化画像情報に対して、前記相関演算部によって算出された相関値、予め設定されたぼけ情報、および観測雑音情報を用いて、時刻ｎ＋１までの情報による時刻ｎ＋１での前記状態量の最適推定値と、時刻ｎまでの情報による時刻ｎ＋１での前記状態量の最適推定値と、前記劣化画像情報を含む観測量の推定誤差と、の関係を規定するための重み係数を算出する重み係数算出部と、
時刻ｎのみの劣化画像情報に対して、前記重み係数算出部によって算出された重み係数を用いて、時刻ｎまたは時刻ｎ＋１までの情報による時刻ｎ＋１での前記状態量の最適推定値を算出する最適推定値算出部と、を有し、
前記劣化画像情報は、注目画素領域と周辺画素領域とからなる処理対象ブロックの画素情報によって構成され、前記時刻ｎは、前記処理対象ブロックの処理の順番（ｎ番目）を表す、
画像復元装置。
【請求項２】
前記相関演算部は、第１の相関演算部および第２の相関演算部を有し、
前記最適推定値算出部は、第１の最適推定値算出部および第２の最適推定値算出部を有し、
前記第１の相関演算部は、
時刻ｎのみの劣化画像情報に対して、時刻ｎまでの情報により前記原画像情報を含む時刻ｎ＋１でのシステムの状態量を推定した場合の推定誤差の相関値行列を算出し、
前記重み係数算出部は、
時刻ｎのみの劣化画像情報に対して、前記第１の相関演算部によって算出された推定誤差の相関値行列、前記ぼけ情報、および前記観測雑音情報を用いて、時刻ｎ＋１までの情報による時刻ｎ＋１での前記状態量の最適推定値と、時刻ｎまでの情報による時刻ｎ＋１での前記状態量の最適推定値と、前記劣化画像情報を含む観測量の推定誤差と、の関係を規定するための重み係数行列を算出し、
前記第１の最適推定値算出部は、
時刻ｎのみの劣化画像情報に対して、時刻ｎまでの情報による時刻ｎ＋１での前記状態量の最適推定値ベクトルを算出し、
前記第２の最適推定値算出部は、
時刻ｎのみの劣化画像情報に対して、前記重み係数算出部によって算出された重み係数行列および前記ぼけ情報を用いて、時刻ｎ＋１までの情報による時刻ｎ＋１での前記状態量の最適推定値ベクトルを算出し、
前記第２の相関演算部は、
時刻ｎのみの劣化画像情報に対して、前記重み係数算出部によって算出された重み係数行列および前記ぼけ情報を用いて、時刻ｎ＋１までの情報により時刻ｎ＋１の前記状態量を推定した場合の推定誤差の相関値行列を算出する、
請求項１記載の画像復元装置。
【請求項３】
前記第１の相関演算部は、
所定の状態遷移行列、与えられた駆動源ベクトルの共分散の要素値、および与えられたまたは前回前記第２の相関演算部によって算出された推定誤差の相関値行列を用いて、前記推定誤差の相関値行列の算出を行い、
前記重み係数算出部は、
前記第１の相関演算部によって算出された推定誤差の相関値行列、前記ぼけ情報を有する与えられた観測遷移行列、および与えられた観測雑音ベクトルの共分散の要素値を用いて、前記重み係数行列の算出を行い、
前記第１の最適推定値算出部は、
前記状態遷移行列、および、与えられたまたは前回前記第２の最適推定値算出部によって算出された状態量の最適推定値ベクトルを用いて、前記状態量の最適推定値ベクトルの算出を行い、
前記第２の最適推定値算出部は、
前記第１の最適推定値算出部によって算出された状態量の最適推定値ベクトル、前記重み係数算出部によって算出された重み係数行列、前記観測遷移行列、および時刻ｎ＋１のみにおける観測量を用いて、前記状態量の最適推定値ベクトルの算出を行い、
前記第２の相関演算部は、
前記重み係数算出部によって算出された重み係数行列、前記観測遷移行列、および前記第１の相関演算部によって算出された推定誤差の相関値行列を用いて、前記推定誤差の相関値行列の算出を行う、
請求項２記載の画像復元装置。
【請求項４】
原画像情報に不必要な情報が混在した劣化画像情報のみから前記原画像情報を推定する画像復元方法であって、
時刻ｎのみの劣化画像情報に対して、時刻ｎまたは時刻ｎ＋１までの情報により前記原画像情報を含む時刻ｎ＋１でのシステムの状態量を推定した場合の推定誤差の相関値を算出する相関演算工程と、
時刻ｎのみの劣化画像情報に対して、前記相関演算工程で算出した相関値、予め設定されたぼけ情報、および観測雑音情報を用いて、時刻ｎ＋１までの情報による時刻ｎ＋１での前記状態量の最適推定値と、時刻ｎまでの情報による時刻ｎ＋１での前記状態量の最適推定値と、前記劣化画像情報を含む観測量の推定誤差と、の関係を規定するための重み係数を算出する重み係数算出工程と、
時刻ｎのみの劣化画像情報に対して、前記重み係数算出工程で算出した重み係数を用いて、時刻ｎまたは時刻ｎ＋１までの情報による時刻ｎ＋１での前記状態量の最適推定値を算出する最適推定値算出工程と、を有し、
前記劣化画像情報は、注目画素領域と周辺画素領域とからなる処理対象ブロックの画素情報によって構成され、前記時刻ｎは、前記処理対象ブロックの処理の順番（ｎ番目）を表す、
画像復元方法。
【請求項５】
前記相関演算工程は、第１の相関演算工程および第２の相関演算工程を有し、
前記最適推定値算出工程は、第１の最適推定値算出工程および第２の最適推定値算出工程を有し、
前記第１の相関演算工程は、
時刻ｎのみの劣化画像情報に対して、時刻ｎまでの情報により前記原画像情報を含む時刻ｎ＋１でのシステムの状態量を推定した場合の推定誤差の相関値行列を算出し、
前記重み係数算出工程は、
時刻ｎのみの劣化画像情報に対して、前記第１の相関演算工程で算出した推定誤差の相関値行列、前記ぼけ情報、および前記観測雑音情報を用いて、時刻ｎ＋１までの情報による時刻ｎ＋１での前記状態量の最適推定値と、時刻ｎまでの情報による時刻ｎ＋１での前記状態量の最適推定値と、前記劣化画像情報を含む観測量の推定誤差と、の関係を規定するための重み係数行列を算出し、
前記第１の最適推定値算出工程は、
時刻ｎのみの劣化画像情報に対して、時刻ｎまでの情報による時刻ｎ＋１での前記状態量の最適推定値ベクトルを算出し、
前記第２の最適推定値算出工程は、
時刻ｎのみの劣化画像情報に対して、前記重み係数算出工程で算出した重み係数行列および前記ぼけ情報を用いて、時刻ｎ＋１までの情報による時刻ｎ＋１での前記状態量の最適推定値ベクトルを算出し、
前記第２の相関演算工程は、
時刻ｎのみの劣化画像情報に対して、前記重み係数算出部によって算出された重み係数行列および前記ぼけ情報を用いて、時刻ｎ＋１までの情報により時刻ｎ＋１の前記状態量を推定した場合の推定誤差の相関値行列を算出する、
請求項４記載の画像復元方法。
【請求項６】
前記第１の相関演算工程は、
所定の状態遷移行列、与えられた駆動源ベクトルの共分散の要素値、および与えられたまたは前回前記第２の相関演算工程で算出した推定誤差の相関値行列を用いて、前記推定誤差の相関値行列の算出を行い、
前記重み係数算出工程は、
前記第１の相関演算工程で算出した推定誤差の相関値行列、前記ぼけ情報を有する与えられた観測遷移行列、および与えられた観測雑音ベクトルの共分散の要素値を用いて、前記重み係数行列の算出を行い、
前記第１の最適推定値算出工程は、
前記状態遷移行列、および、与えられたまたは前回前記第２の最適推定値算出工程で算出した状態量の最適推定値ベクトルを用いて、前記状態量の最適推定値ベクトルの算出を行い、
前記第２の最適推定値算出工程は、
前記第１の最適推定値算出工程で算出した状態量の最適推定値ベクトル、前記重み係数算出工程で算出した重み係数行列、前記観測遷移行列、および時刻ｎ＋１のみにおける観測量を用いて、前記状態量の最適推定値ベクトルの算出を行い、
前記第２の相関演算工程は、
前記重み係数算出工程で算出した重み係数行列、前記観測遷移行列、および前記第１の相関演算工程で算出した推定誤差の相関値行列を用いて、前記推定誤差の相関値行列の算出を行う、
請求項５記載の画像復元方法。
【請求項７】
原画像情報に不必要な情報が混在した劣化画像情報のみから前記原画像情報を推定するための画像復元プログラムであって、
コンピュータに、
時刻ｎのみの劣化画像情報に対して、時刻ｎまたは時刻ｎ＋１までの情報により前記原画像情報を含む時刻ｎ＋１でのシステムの状態量を推定した場合の推定誤差の相関値を算出する相関演算工程と、
時刻ｎのみの劣化画像情報に対して、前記相関演算工程で算出した相関値、予め設定されたぼけ情報、および観測雑音情報を用いて、時刻ｎ＋１までの情報による時刻ｎ＋１での前記状態量の最適推定値と、時刻ｎまでの情報による時刻ｎ＋１での前記状態量の最適推定値と、前記劣化画像情報を含む観測量の推定誤差と、の関係を規定するための重み係数を算出する重み係数算出工程と、
時刻ｎのみの劣化画像情報に対して、前記重み係数算出工程で算出した重み係数を用いて、時刻ｎまたは時刻ｎ＋１までの情報による時刻ｎ＋１での前記状態量の最適推定値を算出する最適推定値算出工程と、
ここで、前記劣化画像情報は、注目画素領域と周辺画素領域とからなる処理対象ブロックの画素情報によって構成され、前記時刻ｎは、前記処理対象ブロックの処理の順番（ｎ番目）を表す、
を実行させるための画像復元プログラム。
【請求項８】
前記相関演算工程は、第１の相関演算工程および第２の相関演算工程を有し、
前記最適推定値算出工程は、第１の最適推定値算出工程および第２の最適推定値算出工程を有し、
前記第１の相関演算工程は、
時刻ｎのみの劣化画像情報に対して、時刻ｎまでの情報により前記原画像情報を含む時刻ｎ＋１でのシステムの状態量を推定した場合の推定誤差の相関値行列を算出し、
前記重み係数算出工程は、
時刻ｎのみの劣化画像情報に対して、前記第１の相関演算工程で算出した推定誤差の相関値行列、前記ぼけ情報、および前記観測雑音情報を用いて、時刻ｎ＋１までの情報による時刻ｎ＋１での前記状態量の最適推定値と、時刻ｎまでの情報による時刻ｎ＋１での前記状態量の最適推定値と、前記劣化画像情報を含む観測量の推定誤差と、の関係を規定するための重み係数行列を算出し、
前記第１の最適推定値算出工程は、
時刻ｎのみの劣化画像情報に対して、時刻ｎまでの情報による時刻ｎ＋１での前記状態量の最適推定値ベクトルを算出し、
前記第２の最適推定値算出工程は、
時刻ｎのみの劣化画像情報に対して、前記重み係数算出工程で算出した重み係数行列および前記ぼけ情報を用いて、時刻ｎ＋１までの情報による時刻ｎ＋１での前記状態量の最適推定値ベクトルを算出し、
前記第２の相関演算工程は、
時刻ｎのみの劣化画像情報に対して、前記重み係数算出部によって算出された重み係数行列および前記ぼけ情報を用いて、時刻ｎ＋１までの情報により時刻ｎ＋１の前記状態量を推定した場合の推定誤差の相関値行列を算出する、
請求項７記載の画像復元プログラム。
【請求項９】
前記第１の相関演算工程は、
所定の状態遷移行列、与えられた駆動源ベクトルの共分散の要素値、および与えられたまたは前回前記第２の相関演算工程で算出した推定誤差の相関値行列を用いて、前記推定誤差の相関値行列の算出を行い、
前記重み係数算出工程は、
前記第１の相関演算工程で算出した推定誤差の相関値行列、前記ぼけ情報を有する与えられた観測遷移行列、および与えられた観測雑音ベクトルの共分散の要素値を用いて、前記重み係数行列の算出を行い、
前記第１の最適推定値算出工程は、
前記状態遷移行列、および、与えられたまたは前回前記第２の最適推定値算出工程で算出した状態量の最適推定値ベクトルを用いて、前記状態量の最適推定値ベクトルの算出を行い、
前記第２の最適推定値算出工程は、
前記第１の最適推定値算出工程で算出した状態量の最適推定値ベクトル、前記重み係数算出工程で算出した重み係数行列、前記観測遷移行列、および時刻ｎ＋１のみにおける観測量を用いて、前記状態量の最適推定値ベクトルの算出を行い、
前記第２の相関演算工程は、
前記重み係数算出工程で算出した重み係数行列、前記観測遷移行列、および前記第１の相関演算工程で算出した推定誤差の相関値行列を用いて、前記推定誤差の相関値行列の算出を行う、
請求項８記載の画像復元プログラム。
【請求項１０】
画像データを入力する画像入力部と、
請求項４から請求項６のいずれかに記載の画像復元方法を実施して、前記画像入力部によって入力された画像データを劣化画像情報として、原画像情報を推定する画像復元処理部と、
前記画像復元処理部によって推定された原画像情報を出力する画像出力部と
を有する画像復元装置。
【請求項１１】
使用者の入力操作により、前記画像復元処理部による画像復元処理に必要な各種設定を行う操作部、
をさらに有する請求項１０記載の画像復元装置。
【請求項１２】
前記操作部は、
前記入力された画像データに対して、前記画像復元処理の対象となる領域を指定する領域指定部、
を有する請求項１１記載の画像復元装置。
【請求項１３】
前記画像復元処理部は、
前記画像復元処理の内容が異なる複数の復元モードを実行可能であり、
前記操作部は、
前記画像復元処理部が実行する復元モードを指定する復元モード指定部、
を有する請求項１１記載の画像復元装置。
【請求項１４】
画像データを入力する画像入力工程と、
前記画像入力工程で入力した画像データを表示する画像表示工程と、
前記画像表示工程で表示した画像データに対して、画像復元処理の対象となる領域を設定する領域設定工程と、
請求項４から請求項６のいずれかに記載の画像復元方法を実施して、前記画像入力工程で入力した画像データのうち、前記領域設定工程で設定した領域に含まれる画像データを劣化画像情報として、原画像情報を推定する画像復元処理工程と、
前記画像復元処理工程で推定した原画像情報を出力する画像出力工程と、
を有する画像復元方法。
【請求項１５】
前記画像復元処理工程は、
前記画像復元処理の内容が異なる複数の復元モードを実行可能であり、
指定された復元モードに対応する内容の画像復元処理を実行する、
を有する請求項１４記載の画像復元方法。
【請求項１６】
請求項４から請求項６のいずれかに記載の画像復元方法以外の他の画像復元方法を実施して、前記画像入力工程で入力した画像データのうち、前記領域設定工程で設定した領域に含まれる画像データを劣化画像情報として、原画像情報を推定する他の画像復元処理工程と、
前記領域設定工程で設定した領域の中から、画像の急激な変化領域を検出する変化領域検出工程と、をさらに有し、
前記変化領域検出工程で検出した画像の急激な変化領域に対して、前記画像復元処理工程を実行し、それ以外の領域に対して、前記他の画像復元処理工程を実行する、
請求項１４記載の画像復元方法。
【発明の詳細な説明】
【技術分野】
【０００１】
本発明は、画像復元装置および画像復元方法に関する。
【背景技術】
【０００２】
近年、画像工学の分野において、劣化した画像から原画像を復元する技術の研究開発が多くなされている。すなわち、原画像（所望情報、クリアな画像）に不必要な情報（ぼけと雑音）が混在した劣化画像（受信情報）から不必要な情報（ぼけと雑音）を取り除き、原画像（所望情報）のみを抽出することは、画像復元の分野において必要不可欠な技術であり、近年盛んに研究開発が行われている。例えば、デジタルカメラ（デジタルスチルカメラおよびデジタルビデオカメラの総称）や携帯電話などで撮影された画像は、手ぶれや焦点ずれなどに起因する「ぼけ」と、暗電流や熱雑音などに起因するガウス性またはインパルス性の「雑音」との影響から、実物と比較して画質が劣化することは避けられない。「画像復元」とは、このように劣化した画像からできるだけ原画像に近い画像を復元することである。【０００３】
現在一般に市場で普及している画像復元技術の多くは、例えば、手ぶれ補正や顔認識、色補正、様々なフィルタなどを用いて、予めぼけと雑音の影響を軽減させる予防的な技術が主流となっている。この結果、近年、特にデジタルカメラの分野では、デジタルカメラの多機能化・高性能化により鮮明な画像が気軽に手に入れられるようになってきている。【０００４】
しかし、これらの予防的な技術は、何度も画像を撮り直せる状況下では問題ないが、古文書などの既に劣化している画像や、スポーツや医療などの分野における瞬間的に変化する画像など、撮り直しが許されない画像に対する復元は、いまだ困難な問題である。ここで、スポーツと医療の分野における瞬間的に変化する画像は、例えば、選手の瞬間的な動作や、肺や心臓などの臓器の瞬間的な様子などである。したがって、今日、特に、撮り直しが許されない状況下における画像復元の重要性が高まってきている。【０００５】
撮り直しが許されない状況下における従来の画像復元技術として、ウィナーフィルタ（Wiener Filter）を用いた画像復元手法が広く知られている（非特許文献１、非特許文献２）。この手法は、フィルタを通して得られる復元画像と原画像との平均２乗誤差を最小にするフィルタであって、最小２乗フィルタとも呼ばれている。この手法は、周波数領域で処理を行うため、確率過程の定常性と半無限長の画像のサイズとを前提にした画像復元手法である。【０００６】
また、他の画像復元技術として、射影フィルタを用いた画像復元手法も知られている（非特許文献３、非特許文献４）。射影フィルタは、原画像と復元画像との近さを直接原画像の空間で評価しており、復元画像の雑音の画像成分が個々の原画像の最良近似、つまり、個々の原画像の正射影となるものの中で、復元画像の雑音成分の２乗平均を最小にするものである。また、この性質から、射影フィルタは、出現頻度に無関係に、最良近似画像を復元する手法である。【０００７】
また、さらに他の画像復元技術として、カルマンフィルタ（Kalman Filter）を用いた画像復元手法も知られている（非特許文献５、非特許文献６）。この手法は、まずステップ１で、劣化画像に対してＡＲ（Auto Regressive：自己回帰）システムのパラメータ（以下「ＡＲ係数」という）を推定した後、ステップ２で、ステップ１で推定したＡＲ係数を用いて状態空間モデル（状態方程式と観測方程式）を構成し、これにカルマンフィルタ理論（カルマンフィルタアルゴリズム）を適用することによって、高性能な画像復元を実現している。
【先行技術文献】
【非特許文献】
【０００８】
【非特許文献１】西宮亮平、外１名、「複数のウィーナーフィルタによる画像復元」（Image Restoration by Using Multiple Wiener Filters）、電子情報通信学会論文誌、Ａ、Ｖｏｌ．Ｊ８３−Ａ、Ｎｏ．７、ｐｐ．８９２−９０２、２０００年７月
【非特許文献２】山根延元、外３名、「画像の局所定常ガウス混合モデルに基づく適応的ウィーナーフィルタを用いた最適雑音除去」（An Optimal Noise Removal Using Adaptive Wiener Filter Based on Locally Stationary Gaussian Mixture Distribution Model for Images）、電子情報通信学会論文誌、Ａ、Ｖｏｌ．Ｊ８５−Ａ、Ｎｏ．９、ｐｐ．９９３−１００４、２００２年９月
【非特許文献３】小川 英光 、外１名、「部分射影フィルタの諸性質」（Properties of Partial Projection Filter）、電子情報通信学会論文誌、Ａ、Ｖｏｌ．Ｊ７１−Ａ、Ｎｏ．２、ｐｐ．５２７−５３４、１９８８年２月
【非特許文献４】小出 裕司、外２名、「信号・画像推定のための射影フィルタ族の統一理論」（A Unified Theory of the Family of Projection Filters for Signal and Image Estimation）、電子情報通信学会論文誌、Ｄ−ＩＩ、Ｖｏｌ．Ｊ７７−Ｄ−ＩＩ、Ｎｏ．７、ｐｐ．１２９３−１３０１、１９９４年７月
【非特許文献５】城崇、外２名、「カルマンフィルタによる高精度画像復元のための画像のモデリングとパラメータ推定」（Image Modeling and Parameter Identification for Image Restoration Using a Kalman Filter）、電子情報通信学会論文誌、Ｄ−ＩＩ、Ｖｏｌ．Ｊ８０−Ｄ−ＩＩ、Ｎｏ．１１、ｐｐ．２９１２−２９１９、１９９７年１１月
【非特許文献６】松村淳、外２名、「効果的な雑音除去のための適応的な画像のモデル化によるカルマンフィルタ」（A Kalman Filter Using Adaptive Image Modeling for Noise Reduction）、電子情報通信学会論文誌、Ｄ−ＩＩ、Ｖｏｌ．Ｊ８６−Ｄ−ＩＩ、Ｎｏ．２、ｐｐ．２１２−２２２、２００３年２月
【発明の概要】
【発明が解決しようとする課題】
【０００９】
しかしながら、ウィナーフィルタを用いた画像復元手法は、画像の劣化状態によらず復元が可能であるという利点を有する反面、非定常性（画像の分散が変化していること）が強い自然画像（加工されていない撮影されたままの画像）に対する復元精度が低いという欠点がある。【００１０】
すなわち、ウィナーフィルタを用いた画像復元手法は、上記のように、周波数領域で処理を行うため、確率過程の定常性と半無限長の画像のサイズとを前提にしているが、現実には、実環境においてこの前提条件が成り立つことは困難であるため、適切に復元されない場合が存在する（自然画像の非定常性の問題）。また、この手法は、最小２乗誤差を評価量とした一括処理であるため、復元画像にぼけが残るという問題が生じる可能性もある（評価量の問題）。【００１１】
具体的には、ウィナーフィルタを用いた画像復元手法では、画像の出現頻度が低くなると、復元精度が低下する、つまり、画像においては、出現頻度が低いエッジ部分の復元に影響するという問題がある。換言すれば、例えば、カメラで撮影された画像は、エッジ部分が多く非定常性が強いため、特にエッジ部分において最適に復元できないという問題がある。ただし、階調や色などが変化していない、つまり、画像の分散が変化していない、単調な背景部分などは、定常性が強いため、通常、良好な復元精度を得ることができる。【００１２】
また、カルマンフィルタを用いた画像復元手法は、ウィナーフィルタを用いた画像復元手法の問題点（自然画像の非定常性の問題と評価量の問題）を解決するための手法であって、ウィナーフィルタを用いた画像復元手法の問題点を解決可能であるという利点を有する反面、処理対象画像（観測画像・劣化画像）にぼけに基づく劣化があると、復元精度が低いという問題がある。【００１３】
すなわち、カルマンフィルタを用いた画像復元手法は、ステップ１で、注目画素とその周囲の画素との相関を考慮すべく、前処理として注目画素とその周囲の画素に対するＡＲ係数を推定した後、ステップ２で、ステップ１で推定したＡＲ係数を用いて状態空間モデル（状態方程式と観測方程式）を構成し、具体的には、ステップ１で推定したＡＲ係数と原画像から状態方程式を構成し、原画像と劣化関数と雑音から観測方程式を構成することにより、カルマンフィルタ理論（カルマンフィルタアルゴリズム）を用いて画像を復元している。したがって、カルマンフィルタを用いた画像復元手法は、定常性を前提としない時間領域のみの処理であり、かつ、推定誤差の分散を評価量とする逐次的処理であるため、ウィナーフィルタを用いた画像復元手法の問題点を解決することができる。【００１４】
しかしながら、一方で、この手法は、ステップ１で推定したＡＲ係数を用いてステップ２のカルマンフィルタアルゴリズムを実行するため、ステップ１でのＡＲ係数の推定精度に劣化画像の復元精度が大きく依存してしまうという問題を有している（ＡＲシステムの問題）。例えば、デジタルカメラの場合、処理対象画像にぼけ（焦点ずれなど）に基づく劣化が含まれると、ステップ１でのＡＲ次数の決定およびＡＲ係数の推定が困難となるため、ステップ２でのカルマンフィルタによる画像復元の精度に影響が出てしまう。【００１５】
この点、一般には、ＡＲ係数の正確な推定は困難である。ＡＲ係数の正確な推定は、例えば、画像復元であれば、クリアな画像（原画像）に依存しているからである。このことは、原画像が既知でなければならないことを意味しているため、リアルタイム処理は困難となる。また、仮に何らかの手法でリアルタイムにＡＲ係数を正確に推定することが可能になったとしても、処理が増加するため演算量の問題は避けられない。また、そもそも、ＡＲ係数の推定はＡＲ係数の次数を決定した後に行われるが、ＡＲ係数の次数を決定するのは非常に困難であり、この点からもＡＲ係数の正確な推定は困難であるといえる。【００１６】
したがって、今日、ウィナーフィルタを用いた画像復元手法の問題点およびカルマンフィルタを用いた画像復元手法の問題点の両方を解決しうる、シンプルで実用的な高性能の画像復元手法、つまり、シンプルな構成を有し、かつ、実環境でも使用可能であり、しかも劣化画像の復元性能（つまり、画像復元能力）が高い画像復元手法が求められている。【００１７】
本発明の目的は、画像復元能力を向上することができるシンプルで実用的な画像復元装置および画像復元方法を提供することである。
【課題を解決するための手段】
【００１８】
本発明の画像復元装置は、原画像情報に不必要な情報が混在した劣化画像情報のみから前記原画像情報を推定する画像復元装置であって、時刻ｎのみの劣化画像情報に対して、時刻ｎまたは時刻ｎ＋１までの情報により前記原画像情報を含む時刻ｎ＋１でのシステムの状態量を推定した場合の推定誤差の相関値を算出する相関演算部と、時刻ｎのみの劣化画像情報に対して、前記相関演算部によって算出された相関値を用いて、時刻ｎ＋１までの情報による時刻ｎ＋１での前記状態量の最適推定値と、時刻ｎまでの情報による時刻ｎ＋１での前記状態量の最適推定値と、前記劣化画像情報を含む観測量の推定誤差と、の関係を規定するための重み係数を算出する重み係数算出部と、時刻ｎのみの劣化画像情報に対して、前記重み係数算出部によって算出された重み係数を用いて、時刻ｎまたは時刻ｎ＋１までの情報による時刻ｎ＋１での前記状態量の最適推定値を算出する最適推定値算出部と、を有する構成を採る。【００１９】
本発明の画像復元装置は、好ましくは、原画像情報に不必要な情報が混在した劣化画像情報のみから前記原画像情報を推定する画像復元装置であって、時刻ｎのみの劣化画像情報に対して、時刻ｎまでの情報により前記原画像情報を含む時刻ｎ＋１でのシステムの状態量を推定した場合の推定誤差の相関値行列を算出する第１の相関演算部と、時刻ｎのみの劣化画像情報に対して、前記第１の相関演算部によって算出された推定誤差の相関値行列を用いて、時刻ｎ＋１までの情報による時刻ｎ＋１での前記状態量の最適推定値と、時刻ｎまでの情報による時刻ｎ＋１での前記状態量の最適推定値と、前記劣化画像情報を含む観測量の推定誤差と、の関係を規定するための重み係数行列を算出する重み係数算出部と、時刻ｎのみの劣化画像情報に対して、時刻ｎまでの情報による時刻ｎ＋１での前記状態量の最適推定値ベクトルを算出する第１の最適推定値算出部と、時刻ｎのみの劣化画像情報に対して、前記重み係数算出部によって算出された重み係数行列を用いて、時刻ｎ＋１までの情報による時刻ｎ＋１での前記状態量の最適推定値ベクトルを算出する第２の最適推定値算出部と、時刻ｎのみの劣化画像情報に対して、時刻ｎ＋１までの情報により時刻ｎ＋１での前記状態量を推定した場合の推定誤差の相関値行列を算出する第２の相関演算部と、を有する構成を採る。【００２０】
本発明の画像復元方法は、原画像情報に不必要な情報が混在した劣化画像情報のみから前記原画像情報を推定する画像復元方法であって、時刻ｎのみの劣化画像情報に対して、時刻ｎまたは時刻ｎ＋１までの情報により前記原画像情報を含む時刻ｎ＋１でのシステムの状態量を推定した場合の推定誤差の相関値を算出する相関演算工程と、時刻ｎのみの劣化画像情報に対して、前記相関演算工程で算出した相関値を用いて、時刻ｎ＋１までの情報による時刻ｎ＋１での前記状態量の最適推定値と、時刻ｎまでの情報による時刻ｎ＋１での前記状態量の最適推定値と、前記劣化画像情報を含む観測量の推定誤差と、の関係を規定するための重み係数を算出する重み係数算出工程と、時刻ｎのみの劣化画像情報に対して、前記重み係数算出工程で算出した重み係数を用いて、時刻ｎまたは時刻ｎ＋１までの情報による時刻ｎ＋１での前記状態量の最適推定値を算出する最適推定値算出工程と、を有するようにした。【００２１】
本発明の画像復元方法は、好ましくは、原画像情報に不必要な情報が混在した劣化画像情報のみから前記原画像情報を推定する画像復元方法であって、時刻ｎのみの劣化画像情報に対して、時刻ｎまでの情報により前記原画像情報を含む時刻ｎ＋１でのシステムの状態量を推定した場合の推定誤差の相関値行列を算出する第１の相関演算工程と、時刻ｎのみの劣化画像情報に対して、前記第１の相関演算工程で算出した推定誤差の相関値行列を用いて、時刻ｎ＋１までの情報による時刻ｎ＋１での前記状態量の最適推定値と、時刻ｎまでの情報による時刻ｎ＋１での前記状態量の最適推定値と、前記劣化画像情報を含む観測量の推定誤差と、の関係を規定するための重み係数行列を算出する重み係数算出工程と、時刻ｎのみの劣化画像情報に対して、時刻ｎまでの情報による時刻ｎ＋１での前記状態量の最適推定値ベクトルを算出する第１の最適推定値算出工程と、時刻ｎのみの劣化画像情報に対して、前記重み係数算出工程で算出した重み係数行列を用いて、時刻ｎ＋１までの情報による時刻ｎ＋１での前記状態量の最適推定値ベクトルを算出する第２の最適推定値算出工程と、時刻ｎのみの劣化画像情報に対して、時刻ｎ＋１までの情報により時刻ｎ＋１での前記状態量を推定した場合の推定誤差の相関値行列を算出する第２の相関演算工程と、を有するようにした。
【発明の効果】
【００２２】
本発明によれば、画像復元能力を向上することができるシンプルで実用的な画像復元装置および画像復元方法を得ることができる。
【図面の簡単な説明】
【００２３】
【図１】本発明の一実施の形態に係る画像復元装置の構成を示すブロック図
【図２】図２Ａは、図１の画像復元処理部の構成を示すブロック図、図２Ｂは、図２Ａの第１復元処理部の構成を示すブロック図
【図３】画像の劣化の過程を説明するための図
【図４】発明手法の状態空間モデルのシステム構成を表すブロック線図
【図５】発明手法の状態方程式の定式化の具体例を説明するための図であり、特に、図５Ａは、状態空間モデルの状態過程を示す図、図５Ｂは、処理対象ブロックおよびその時間変化の一例を示す図、図５Ｃは、状態方程式の具体例を示す図
【図６】発明手法の状態方程式の構成を説明するための図
【図７】発明手法の観測方程式の定式化の具体例を説明するための図であり、特に、図７Ａは、状態空間モデルの観測過程を示す図、図７Ｂは、観測方程式の具体例を示す図
【図８】従来の一般的な観測方程式の構成（一般例）を説明するための図であり、特に、図８Ａは、従来の観測方程式の定義を示す図、図８Ｂは、観測量に影響する状態量の範囲を視覚的に示す図
【図９】発明手法の観測方程式の構成を説明するための図であり、特に、図９Ａは、発明手法の観測方程式の定義を示す図、図９Ｂは、観測量に影響する状態量の範囲を視覚的に示す図
【図１０】観測遷移行列を構成する係数の割り当て方法を説明するための図
【図１１】図９Ａに示す観測遷移行列に基づく、各観測量に影響する状態量の範囲を、割り当てられた係数と共に示す図
【図１２】発明手法のアルゴリズムの一例を示す図
【図１３】図１２のアルゴリズムを実行する処理手順を示すフローチャート
【図１４】発明手法を視覚的にまとめた説明図
【図１５】シミュレーション条件を説明するための図
【図１６】原画像「カメラマン」を説明するための図
【図１７】原画像「カメラマン」に対するシミュレーション結果（視覚評価）を示す図
【図１８】図１７中の破線の円の周辺領域を拡大した図
【図１９】原画像「カメラマン」に対するシミュレーション結果（視覚評価）を示す図
【図２０】図１９中の破線の円の周辺領域を拡大した図
【図２１】原画像「レナ」に対するシミュレーション結果（視覚評価）を示す図
【図２２】図２１中の破線の円の周辺領域を拡大した図
【図２３】原画像に対するシミュレーション結果（客観評価）を示す図
【図２４】原画像に対するシミュレーション結果（主観評価）を示す図
【図２５】図１の画像復元装置の動作の一例を示すフローチャート
【図２６】図２５のステップＳ２３００の検出方法の一例を説明するための図
【図２７】復元モードの実装の一例を説明するための図
【図２８】復元モードの実装の他の例を説明するための図
【発明を実施するための形態】
【００２４】
以下、本発明の実施の形態について、図面を参照して詳細に説明する。
【００２５】
図１は、本発明の一実施の形態に係る画像復元装置の構成を示すブロック図である。ここでは、本発明の画像復元方法を適用した好適な画像復元装置について、特に各種の用途に対応可能な汎用的な画像復元装置を例にとって説明する。
【００２６】
図１に示す画像復元装置１００は、コンピュータで構成されており、大別して、画像入力装置１１０、入力インタフェース部１２０、操作部１３０、内部インタフェース部１４０、記憶部１５０、画像復元処理部１６０、出力インタフェース部１７０、および画像出力装置１８０を有する。
【００２７】
画像入力装置１１０は、復元処理の対象となる画像データ（劣化画像）をデジタルデータとしてコンピュータに入力するための入力装置である。入力する画像は、静止画像でも動画像でもよい。画像入力装置１１０としては、例えば、カメラ１１２やスキャナ１１４、記録メディア１１６、モデム１１８などを用いることができる。カメラ１１２は、撮像機能を有するすべての装置を意味し、例えば、デジタルカメラ（デジタルスチルカメラとデジタルビデオカメラ）のほかに、カメラ機能を搭載した携帯電話や、防犯カメラ（監視カメラ）、画像診断を行うための医療機器（内視鏡、レントゲン、エコー、ＣＴ、ＭＲＩなど）などを含むことができる。スキャナ１１４は、代表的な画像入力装置の１つであり、ネガやポジといったファイルからの読み込みを専用に行うフィルムスキャナなども含んでいる。記録メディア１１６は、画像データを記録可能な記録媒体を広く意味し、例えば、磁気ディスク（ＨＤＤやＦＤなど）や、光ディスク（ＣＤやＤＶＤ、ＢＤなど）、光磁気ディスク（ＭＯ）、フラッシュメモリ（メモリカードやＵＳＢメモリなど）などである。モデム１１８は、外部の通信ネットワーク（例えば、電話回線やＬＡＮ、インターネットなど）と接続するための装置である。なお、画像入力装置１１０の種類は、画像復元装置１００の用途に応じて適宜選択すればよい。
【００２８】
入力インタフェース部１２０は、画像入力装置１１０から提供される画像データをコンピュータで処理可能なデータ形式に変換するなどの入力処理を行う。図示しないが、入力インタフェース部１２０は、画像入力装置１１０の種類に応じて別個独立に設けられる。例えば、記録メディア１１６の入力インタフェース部１２０は、ドライブとも呼ばれ、記録メディア１１６の種類に応じて様々な種類のドライブが使用可能である。なお、ドライブは、記録メディアを読み書きする装置であり、記録メディアに関する限り、通常は、入力インタフェース部１２０と出力インタフェース部１７０は一体化されている。また、通常、モデム１１８は、画像入力装置１１０としても画像出力装置１８０としても機能しうるため、モデム１１８に関しても、通常、入力インタフェース部１２０と出力インタフェース部１７０は一体化されている。入力インタフェース部１２０は、コンピュータ本体の内部に格納された内蔵カード（ボード）であってもよいし、内部インタフェース部１４０を介して接続された外部設置型機器であってもよい。
【００２９】
なお、画像入力装置１１０が画像情報をアナログデータとして出力する場合、対応する入力インタフェース部１２０は、サンプリング部およびＡ／Ｄ変換部（共に図示せず）を有する。サンプリング部は、所定のサンプリング周波数で、入力されたアナログ信号をサンプリング処理し、Ａ／Ｄ変換部に出力する。サンプリング周波数は、復元処理対象（情報源）の種類に応じて変更可能である。Ａ／Ｄ変換部は、サンプリングされた信号の振幅値を所定の分解能でＡ／Ｄ変換処理する。
【００３０】
操作部１３０は、典型的にはキーボートやマウス、タッチパネルなどであるが、音声認識装置などを用いてもよい。使用者は、操作部１３０を用い、例えば、後述するディスプレイ１８２で確認をしながらコンピュータを操作することができる。また、操作部１３０は、パラメータ設定部１３２、領域指定部１３４、および復元モード指定部１３６を有する。パラメータ設定部１３２は、使用者の入力操作により、後で詳述する本実施の形態における画像復元処理に必要な各種パラメータの値を設定する。領域指定部１３４は、使用者の入力操作により、入力された画像に対して画像復元処理の対象となる領域（画像の特定の範囲）を指定する。復元モード指定部１３６は、使用者の入力操作により、後述する復元モードを指定する。
【００３１】
内部インタフェース部１４０は、コンピュータ本体の内部にあって、入力インタフェース部１２０、操作部１３０、記憶部１５０、画像復元処理部１６０、および出力インタフェース部１７０を相互に接続する機能を有する。コンピュータの内部における各種信号のやり取りは、内部インタフェース部１４０を介して行われる。
【００３２】
記憶部１５０は、主記憶装置１５２および補助記憶装置１５４を有する。主記憶装置１５２は、コンピュータの本体を構成する一要素であって、主にプログラムとデータを格納する。補助記憶装置１５４は、主記憶装置１５２の容量の不足を補う記憶装置である。補助記憶装置１５４は、典型的にはハードディスク（ＨＤ）であるが、ＣＤ−ＲＯＭやＤＶＤ、ＳＳＤ（Solid State Drive）、フラッシュメモリなどの可搬性のあるものであってもよいし、また、それらの組み合わせであってもよい。本実施の形態における画像復元処理を実行させるプログラム（画像復元アルゴリズム）は、予め記憶部１５０（主記憶装置１５２または補助記憶装置１５４）に格納されていてもよいし、記録メディア１１６から両インタフェース部１２０、１４０を介して記憶部１５０にインストールしたり、モデム１１８および両インタフェース部１２０、１４０を介して外部から記憶部１５０にダウンロードしたりしてもよい。
【００３３】
なお、画像入力装置１１０から画像データを取り込み、取り込んだ画像データに対して画像復元処理を行い、復元処理後の画像データを画像出力装置１８０から取り出すという一連の処理を実行するためには、取り込んだ画像データを格納しておく記憶領域と、データ処理の上で一時的に必要な記憶領域（作業領域または作業用メモリと呼ばれる）と、出力すべき画像データを格納しておく記憶領域とが必要である。これらの記憶領域は、主記憶装置１５２または補助記憶装置１５４の上に取られるが、ここでは、説明の便宜上、復元処理後の画像データを後述するディスプレイ１８２に出力する場合を想定して、特に表示メモリ１５６を図面上に別途表記している。
【００３４】
画像復元処理部１６０は、本発明の特徴的な構成要素であり、後述する内蔵された画像復元アルゴリズムを実行する。カルマンフィルタを用いた従来の画像復元手法では、２段階処理、つまり、ＡＲ次数を決定し、ＡＲ係数を推定した後、この推定したＡＲ係数を用いて状態空間モデル（状態方程式と観測方程式）を構成してカルマンフィルタを実行することによって、画像復元を実現しているのに対し、本発明の画像復元方法（以下「発明手法」という）では、状態方程式と観測方程式から構成される新しい予測法によって画像復元を実現している。具体的には、発明手法では、ＡＲシステムのコンセプトを必要としない新しい状態空間モデル（状態方程式と観測方程式）を構成することによって、より具体的には、状態方程式をクリアな画像情報（原画像情報）のみで構成し、かつ、観測方程式を劣化画像情報、クリアな画像情報（原画像情報）、ぼけ情報、および雑音で構成した新たな状態空間モデルを用いることによって、画像復元を実現している。画像復元の対象となる画像は、静止画像でも動画像でもよい。なお、本発明の画像復元方法については、後で詳細に説明する。
【００３５】
画像出力装置１８０は、コンピュータ（画像復元処理部１６０）で復元処理された画像データ（復元画像）を所定の出力形態で外部に取り出すための出力装置である。画像出力装置１８０としては、例えば、ディスプレイ１８２やプリンタ１８４、記録メディア１８６、モデム１８８などを用いることができる。記録メディア１８６およびモデム１８８は、画像入力装置１１０としての記録メディア１１６およびモデム１１８とそれぞれ共用してもよい。なお、画像出力装置１８０の種類は、画像復元装置１００の用途に応じて適宜選択すればよい。
【００３６】
出力インタフェース部１７０は、コンピュータ（画像復元処理部１６０）で復元処理された画像データ（復元画像）を画像出力装置１８０に出力可能なデータ形式に変換するなどの出力処理を行う。図示しないが、出力インタフェース部１７０は、画像出力装置１８０の種類に応じて別個独立に設けられる。上記のように、記録メディアおよびモデムに関しては、通常、入力インタフェース部１２０と出力インタフェース部１７０は一体化されている。出力インタフェース部１７０も、入力インタフェース部１２０と同様に、コンピュータ本体の内部に格納された内蔵カード（ボード）であってもよいし、内部インタフェース部１４０を介して接続された外部設置型機器であってもよい。
【００３７】
図２Ａは、図１の画像復元処理部１６０の構成を示すブロック図であり、図２Ｂは、図２Ａの第１復元処理部１６０ａの構成を示すブロック図である。
【００３８】
本実施の形態では、図２Ａに示すように、画像復元処理部１６０は、第１復元処理部１６０ａと第２復元処理部１６０ｂを有する。第１復元処理部１６０ａは、本発明の画像復元方法を実施する。第２復元処理部１６０ｂは、本発明の画像復元方法以外の画像復元方法、例えば、ウィナーフィルタを用いた画像復元手法やカルマンフィルタを用いた画像復元手法などの従来の画像復元方法を実施する。これによって、第１復元処理部１６０ａと第２復元処理部１６０ｂとを用いた画像復元処理が可能となる。なお、以下では、第２復元処理部１６０ｂで実施する従来の画像復元方法として、例えば、ウィナーフィルタを用いた画像復元手法（以下「従来手法」という）を例にとる。
【００３９】
第１復元処理部１６０ａは、図２Ｂに示すように、初期設定部１６２、相関演算部１６４、重み係数算出部１６６、および最適推定値算出部１６８を有する。第１復元処理部１６０ａは、各部１６２〜１６８の協働の下、内蔵された画像復元アルゴリズム（発明手法）を実行することによって、取り込んだ画像データ（劣化画像）のみからクリアな画像（原画像）を復元する。このとき、初期設定部１６２は、発明手法の画像復元アルゴリズムの初期設定を行い、相関演算部１６４は、原画像（所望情報、クリアな画像）の推定誤差の相関演算を行い、重み係数算出部１６６は、原画像（所望情報）の最適推定値の算出に必要な重み係数の算出を行い、最適推定値算出部１６８は、原画像（所望情報）の最適推定値の算出を行う。なお、各部１６２〜１６８の具体的な処理内容は、後で詳細に説明する。
【００４０】
以下、第１復元処理部１６０ａで行われる画像復元処理動作について詳細に説明するが、ここでは、まず、画像復元の理論となる劣化画像モデルについて説明しておく。
【００４１】
劣化した画像は、一般的に、原画像とぼけの点拡がり関数（ＰＳＦ：Point Spread Function）の畳み込み（コンボリューション）に雑音を加えたモデルによって得られる。すなわち、イメージングシステム（対象物を撮像する装置や画像を生成するシステムなど）によって得られる劣化画像ｇ（ｘ，ｙ）は、対象物を表す原画像をｆ（ｎ，ｍ）、ぼけの点拡がり関数をｈ（ｘ，ｙ）、雑音をｎ（ｘ，ｙ）でそれぞれ表すと、次の式（１）で表される。なお、ぼけの点拡がり関数ｈ（ｘ，ｙ）は、撮影条件なども含めたイメージングシステムの特性を表している。
【００４２】
【数１】
0005447867.tif 000002
【００４３】
図３は、画像の劣化の過程を説明するための図である。
【００４４】
例えば、図３の左端に示す原画像にぼけが生じると、図３の中央に示す画像となり、この画像にさらに雑音が加わると、図３の右端に示す劣化画像となる。ぼけは、ある画素がその周囲の画素に影響することによって生じるが、雑音は、画素に関係なく発生する。上記のように、カメラなどで撮影された画像の場合、ぼけは、手ぶれや焦点ずれなどに起因して生じ、雑音は、暗電流や熱雑音などに起因して不可避的に発生する。図３に示すように、劣化画像の復元は、劣化画像からぼけと雑音を取り除くことにほかならない。
【００４５】
発明手法では、上記のように、ＡＲシステムのコンセプトを用いないようにするために、新しい状態空間モデル（状態方程式と観測方程式）を構成する。すなわち、クリアな画像情報（原画像情報）のみを用いて状態方程式を構成し、かつ、劣化画像情報、クリアな画像情報（原画像情報）、ぼけ情報、および雑音を用いて観測方程式を構成する。
【００４６】
具体的には、発明手法では、新しい状態空間モデル（状態方程式と観測方程式）を構成し、この新しい状態空間モデルは、次の式（２）で表される。ただし、状態方程式において、ベクトルｘｐ１（ｎ）は状態ベクトル（原画像情報）、行列Φｐ１は状態遷移行列、ベクトルδｐ１（ｎ）は駆動源である。また、観測方程式において、ベクトルｙｐ１（ｎ）は観測ベクトル（観測画像情報、劣化画像情報）、行列Ｍｐ１は観測遷移行列、ベクトルεｐ１（ｎ）は観測雑音である。また、「ｎ」とは、装置の時刻ｎである。時刻ｎは、注目画素を含む複数の周辺画素からなる処理対象ブロックの処理の順番（ｎ番目）を表している。なお、ここでは、発明手法として１つの具体的な手法（以下「発明手法１」という）を提示する。以下の説明において、添え字「ｐ１」は、発明手法１に係るものであることを示している。
【００４７】
【数２】
0005447867.tif 000003
【００４８】
図４は、この状態空間モデルにおけるシステム構成をブロック線図で表したものである。
【００４９】
図４に示すように、この状態空間モデルは、状態過程と観測過程とで構成される。状態過程は、状態方程式によって記述され、観測過程は、観測方程式によって記述される。図４において、「２０１」は時刻ｎにおける状態ベクトルｘｐ１（ｎ）、「２０２」は時刻ｎ＋１における状態ベクトルｘｐ１（ｎ＋１）、「２０３」は時刻ｎにおける観測ベクトルｙｐ１（ｎ）、「２０４」は時刻ｎにおける観測雑音ベクトルεｐ１（ｎ）、「２０５」は時刻ｎ＋１における駆動源ベクトルδｐ１（ｎ＋１）、「２０６」は状態遷移行列Φｐ１、「２０７」は観測遷移行列Ｍｐ１である。式（２）における状態方程式は、観測対象のシステムを状態空間モデルで記述したものであり、内部状態つまり状態変数（ここでは、状態ベクトルｘｐ１（ｎ）、ｘｐ１（ｎ＋１））の時間に対する生成過程を表している。また、式（２）における観測方程式は、何らかの観測装置を通じて観測する過程を記述したものであり、観測結果（ここでは、観測ベクトルｙｐ１（ｎ））が、被観測量つまり入力（ここでは、状態ベクトルｘｐ１（ｎ））に依存して生成される様子を示している。
【００５０】
図５は、発明手法の状態方程式の定式化の具体例を説明するための図であり、特に、図５Ａは、状態空間モデルの状態過程を示す図、図５Ｂは、処理対象ブロックおよびその時間変化の一例を示す図、図５Ｃは、状態方程式の具体例を示す図である。
【００５１】
本実施の形態では、劣化画像の復元処理は、処理対象画素のみを用いるのではなく、その周辺の画素も含めて処理を行う。すなわち、Ｊ×Ｊの大きさの画像のある処理対象画素を中心にして、その周辺Ｋ×Ｋ（ただし、Ｊ＞Ｋ）の領域（以下「注目領域」という）を処理対象ブロックとし、この処理対象ブロックの中の画素情報を用いて、画像の復元処理を行う。このように、注目領域とは、画像復元においてＫ×Ｋの画素を用いて処理を行う範囲を意味する。
【００５２】
例えば、Ｊ＝２５６、Ｋ＝３とした場合、図５に示すように、２５６×２５６の大きさの画像のある処理対象画素を中心にして、その周辺３×３の注目領域が処理対象ブロックとなる。図面上、処理対象ブロックは、網かけされている。画像復元処理は、この３×３の注目領域の中の画素情報を用いて行われる。単に説明の便宜上、図５Ｂに示すように、２５６×２５６の大きさの画像の一部の画素に１〜３６の番号を付したとすると、時刻ｎにおける状態ベクトルｘｐ１（ｎ）に対応する処理対象ブロックは、番号が１、２、３、７、８、９、１３、１４、１５の３×３の９個の画素を構成要素とし、次の時刻ｎ＋１における状態ベクトルｘｐ１（ｎ＋１）に対応する処理対象ブロックは、番号が７、８、９、１３、１４、１５、１９、２０、２１の３×３の９個の画素を構成要素としている。
【００５３】
このとき、式（２）における状態方程式は、図５Ｃにも示す次の式（３）で定義される。ここで、状態ベクトルｘｐ１（ｎ）は、原画像情報を含む状態量として３×３の処理対象ブロックに含まれる９個の画素情報ｘ１（ｎ）、ｘ２（ｎ）、ｘ３（ｎ）、ｘ４（ｎ）、ｘ５（ｎ）、ｘ６（ｎ）、ｘ７（ｎ）、ｘ８（ｎ）、ｘ９（ｎ）を要素とする９×１次のベクトルである。また、状態遷移行列Φｐ１は、式（３）で定義された９×９次の行列であり、駆動源ベクトルδｐ１（ｎ）は、同じく式（３）で定義された９×１次のベクトルである。
【００５４】
【数３】
0005447867.tif 000004
【００５５】
図６は、式（３）で表される状態方程式の構成を説明するための図である。
【００５６】
式（３）で表される状態方程式の特徴は、状態遷移行列Φｐ１の一部の要素を「１」に設定し、残りの要素をすべて「０」に設定した点、および、駆動源ベクトルδｐ１（ｎ）の一部の要素を有色信号である状態量ｘｉ（ｎ）（ｉ＝７，８，９）で表した点にある。これは、図６に示すように、現在の時刻ｎ＋１における処理対象ブロックと過去の時刻ｎにおける処理対象ブロックとが重なる部分において、現在の状態量（原画像の画素情報）ｘ１（ｎ＋１）、ｘ２（ｎ＋１）、ｘ３（ｎ＋１）、ｘ４（ｎ＋１）、ｘ５（ｎ＋１）、ｘ６（ｎ＋１）と過去の状態量（原画像の画素情報）ｘ４（ｎ）、ｘ５（ｎ）、ｘ６（ｎ）、ｘ７（ｎ）、ｘ８（ｎ）、ｘ９（ｎ）とを対応させるためである。この結果、ｘｐ１（ｎ）とｘｐ１（ｎ＋１）の関係を表す状態方程式は、クリアな画像からなるｘｐ１（ｎ）と、０と１からなるΦｐ１と、有色信号であるクリアな画像からなる駆動源ベクトルδｐ１（ｎ）とから構成されるため、所望の状態量（原画像の画素情報）、つまり、クリアな画像情報（原画像情報）のみで構成されることになる。
【００５７】
図７は、発明手法の観測方程式の定式化の具体例を説明するための図であり、特に、図７Ａは、状態空間モデルの観測過程を示す図、図７Ｂは、観測方程式の具体例を示す図である。
【００５８】
図５および図６の例に対応して、式（２）における観測方程式は、図７Ｂにも示す次の式（４）で定義される。ここで、観測ベクトルｙｐ１（ｎ）は、劣化画像情報を含む観測量として３×３の処理対象ブロックに含まれる９個の画素情報ｙ１（ｎ）、ｙ２（ｎ）、ｙ３（ｎ）、ｙ４（ｎ）、ｙ５（ｎ）、ｙ６（ｎ）、ｙ７（ｎ）、ｙ８（ｎ）、ｙ９（ｎ）を要素とする９×１次のベクトルである。また、観測遷移行列Ｍｐ１は、式（４）で定義された９×９次の行列であって、画像の劣化モデルにおけるぼけの点拡がり関数（ＰＳＦ）に対応している。観測遷移行列Ｍｐ１を構成する各要素ｈｒ，ｓ（ｒ，ｓはｈの座標であり、ｒ，ｓ＝−１，０，１）は、既知であって、予めデータ化して適切に設定されている。また、観測雑音ベクトルεｐ１（ｎ）は、３×３の処理対象ブロックに含まれる９個の画素に対応する観測雑音ε１（ｎ）、ε２（ｎ）、ε３（ｎ）、ε４（ｎ）、ε５（ｎ）、ε６（ｎ）、ε７（ｎ）、ε８（ｎ）、ε９（ｎ）を要素とする９×１次のベクトルである。
【００５９】
【数４】
0005447867.tif 000005
【００６０】
図８は、従来の一般的な観測方程式の構成（一般例）を説明するための図であり、図９は、式（４）で表される観測方程式の構成を説明するための図である。特に、図８Ａは、従来一般的な観測方程式の定義を示す図、図８Ｂは、観測量に影響する状態量の範囲を視覚的に示す図である。また、図９Ａは、発明手法の観測方程式の定義を示す図、図９Ｂは、観測量に影響する状態量の範囲を視覚的に示す図である。なお、図８および図９の説明においては、便宜上、観測雑音ベクトルεｐ１（ｎ）は省略する。
【００６１】
画像においてぼけに基づく劣化とは、上記のように、ある画素がその周囲の画素の影響を受けて劣化することである。そこで、従来は、一般的に、図８Ａに示すように観測方程式を定義していた。この場合、ぼけは、周囲のすべての画素の影響を受けることに起因して生じる。すなわち、処理対象ブロックの中の各画素は、処理対象ブロックの中のすべての画素の影響を受けることになる。言い換えれば、各観測量ｙｉ（ｉ＝１，２，…，９）は、すべての状態量ｘｉ（ｉ＝１，２，…，９）の影響を受けて成り立つ。これを視覚的に示すと、例えば、図８Ｂに示すようになる。すなわち、図８Ｂに示すように、例えば、観測量がｙ１の場合もｙ２の場合も、３×３（＝９）個のすべての状態量ｘ１〜ｘ９の影響を受けることになる。
【００６２】
これに対し、式（４）で表される観測方程式の特徴は、観測遷移行列Ｍｐ１の一部の要素を、３×３（＝９）個の係数ｈｒ，ｓ（ｒ，ｓ＝−１，０，１）を用いて規則的に設定し、残りの要素をすべて「０」に設定した点にある。例えば、図９Ａは、観測方程式を展開して得られる一部の式（ｙ１（ｎ）、ｙ２（ｎ）に関する式）を示している。観測遷移行列Ｍｐ１に対する係数ｈｒ，ｓの割り当て方法は、次の通りである。すなわち、９個の係数ｈｒ，ｓ（ｒ，ｓ＝−１，０，１）を図１０に示すようにマトリックス状に配列し（以下「係数マトリックス」という）、係数マトリックスの中央（つまり、係数ｈ０，０）の位置を観測量ｙｉ（ｉ＝１，２，…，９）の注目画素ｉの位置に合わせたときの、その係数マトリックスに従って、観測遷移行列Ｍｐ１の一部の要素に係数ｈｒ，ｓ（ｒ，ｓ＝−１，０，１）を割り当てる。これを視覚的に示すと、例えば、図９Ｂに示すようになる。すなわち、図９Ｂに示すように、例えば、観測量がｙ１の場合は、３×３（＝９）個の状態量ｘ１〜ｘ９のうち、５個のｘ１〜ｘ５の影響を受け、観測量がｙ２の場合は、６個のｘ１〜ｘ６の影響を受けることになる。
【００６３】
図１１は、図９Ａに示す観測遷移行列Ｍｐ１に基づく、各観測量ｙｉに影響する状態量ｘｉの範囲を、割り当てられた係数ｈｒ，ｓ（ｒ，ｓ＝−１，０，１）と共に示す図である。なお、ここでは、説明の便宜上、図１１Ａ〜図１１Ｉの各図において、左側の処理対象ブロックには、注目する観測量ｙｉのみを示し、右側の処理対象ブロックには、影響を与える各状態量ｘｉの場所に、この状態量ｘｉに乗算されるｈｒ，ｓのみを示している。
【００６４】
式（４）で表される観測方程式では、観測量ｙ１は、図１１Ａに示すように、５個の状態量ｘ１〜ｘ５の影響を受け、観測量ｙ２は、図１１Ｂに示すように、６個の状態量ｘ１〜ｘ６の影響を受け、観測量ｙ３は、図１１Ｃに示すように、７個の状態量ｘ１〜ｘ７の影響を受ける。また、観測量ｙ４は、図１１Ｄに示すように、８個の状態量ｘ１〜ｘ８の影響を受け、観測量ｙ５は、図１１Ｅに示すように、９個すべての状態量ｘ１〜ｘ９の影響を受け、観測量ｙ６は、図１１Ｆに示すように、８個の状態量ｘ２〜ｘ９の影響を受ける。また、観測量ｙ７は、図１１Ｇに示すように、７個の状態量ｘ３〜ｘ９の影響を受け、観測量ｙ８は、図１１Ｈに示すように、６個の状態量ｘ４〜ｘ９の影響を受け、観測量ｙ９は、図１１Ｉに示すように、５個の状態量ｘ５〜ｘ９の影響を受ける。
【００６５】
図１２は、発明手法のアルゴリズムの一例を示す図である。なお、このアルゴリズムは、画像の種類を問わず、静止画像のみならず動画像にも適用可能である。
【００６６】
発明手法１のアルゴリズムは、図１２に示すように、初期設定（Initialization）の過程と反復（Iteration）の過程とに大別され、反復過程は、新しい状態空間モデル（状態方程式と観測方程式）を基に構成されている。反復過程では、１〜６の手順を逐次繰り返す。
【００６７】
図１３は、図１２のアルゴリズムを実行する処理手順を示すフローチャートである。
【００６８】
まず、初期設定部１６２で、初期設定を行う（Ｓ１０００）。具体的には、初期設定部１６２において、状態ベクトル、つまり、状態量である所望信号（原画像信号）ベクトルの最適推定値（以下「所望信号の最適推定値ベクトル」という）の初期値ｘｐ１（０｜０）、所望信号ベクトルの推定誤差（以下「所望信号の推定誤差ベクトル」という）の相関行列の初期値Ｐｐ１（０｜０）、駆動源ベクトルの共分散の値Ｒδｐ１（ｎ）［ｉ，ｊ］、および、観測雑音ベクトルの共分散の値Ｒεｐ１（ｎ）［ｉ，ｊ］を、次の式（５）に示すようにそれぞれ設定する。また、図１２では図示しないが、状態遷移行列Φｐ１および観測遷移行列Ｍｐ１を、式（３）および式（４）に示すようにそれぞれ設定する。なお、ここでは、時刻ｎのカウンタの初期値を「０」に設定している。また、ベクトルと行列の要素を示す場合、ベクトルａ（ｎ）のｉ番目の要素をａ（ｎ）［ｉ］と表記し、行列Ａ（ｎ）のｉ行ｊ列の要素をＡ（ｎ）［ｉ，ｊ］と表記することとする。また、図１２において、所望信号の最適推定値ベクトルｘｐ１は、ハットｘｐ１と表記する。
【００６９】
【数５】
0005447867.tif 000006
【００７０】
ただし、ベクトル０Ｋは、Ｋ次元の零ベクトルであり、行列ＩＫは、Ｋ次の単位行列である。また、Ｋは、Ｋ×Ｋの大きさの注目領域のその次元である。Ｋ＝３の場合、Ｋ２−Ｋ＝６となり、「ｉ，ｊ＞Ｋ２−Ｋ」は、「ｉ，ｊ＝７，８，９」に対応している。また、Ｅ［δｐ１（ｎ），δｐ１Ｔ（ｎ）］は、駆動源ベクトルδｐ１（ｎ）の自己相関行列である。Ｅ［εｐ１（ｎ）εｐ１Ｔ（ｎ）］［ｉ，ｊ］は、観測雑音ベクトルεｐ１の自己相関行列であり、ここではσ２ε［ｉ，ｊ］と等しく、既知と仮定している。ここでいう「既知」とは、別の任意の方法（アルゴリズム）で求められて与えられることを意味する。もし観測雑音ベクトルεｐ１（ｎ）が白色雑音であり零平均であれば、σ２εは所定のサンプル数の分散で与えられる。
【００７１】
次に、相関演算部１６４で、ｎ→（ｎ＋１）の誤差共分散行列、つまり、ｎ→（ｎ＋１）の推定誤差の相関値（行列）（以下「相関行列」ともいう）を計算する（Ｓ１１００）。具体的には、相関演算部１６４において、時刻ｎまでの情報により時刻ｎ＋１での所望信号ベクトルを推定した場合の誤差（所望信号の推定誤差ベクトル）の相関行列Ｐｐ１（ｎ＋１｜ｎ）を計算する。この計算は、ステップＳ１０００でそれぞれ設定した状態遷移行列Φｐ１および駆動源ベクトルの共分散Ｒδｐ１（ｎ＋１）［ｉ，ｊ］の値、ならびに、ステップＳ１０００で設定した（ｎ＝０の場合）または前回のステップＳ１６００で計算した（ｎ≧１の場合）所望信号の推定誤差ベクトルの相関行列Ｐｐ１（ｎ｜ｎ）を用いて、次の式（６）により行う。なお、このステップＳ１１００は、図１２の反復過程の手順１に相当する。
【００７２】
【数６】
0005447867.tif 000007
【００７３】
次に、重み係数算出部１６６で、重み係数（行列）の計算を行う（Ｓ１２００）。具体的には、重み係数算出部１６６において、観測量である観測信号ベクトルの推定誤差（以下「観測信号の推定誤差ベクトル」という）に重み係数（行列）をかけて、時刻ｎまでの情報による時刻ｎ＋１での所望信号の最適推定値ベクトルｘｐ１（ｎ＋１｜ｎ）を加えたものが、時刻ｎ＋１までの情報による時刻ｎ＋１での所望信号の最適推定値ベクトルｘｐ１（ｎ＋１｜ｎ＋１）になるような、その重み係数行列Ｋｐ１（ｎ＋１）を計算する。この計算は、ステップＳ１１００で計算した所望信号の推定誤差ベクトルの相関行列Ｐｐ１（ｎ＋１｜ｎ）、ならびに、ステップＳ１０００でそれぞれ設定した観測遷移行列Ｍｐ１および観測雑音ベクトルの共分散Ｒεｐ１（ｎ＋１）［ｉ，ｊ］の値を用いて、次の式（７）により行う。なお、このステップＳ１２００は、図１２の反復過程の手順２に相当する。
【００７４】
【数７】
0005447867.tif 000008
【００７５】
次に、最適推定値算出部１６８で、ｎ→（ｎ＋１）の状態量（所望信号）の最適推定値（ベクトル）を計算する（Ｓ１３００）。具体的には、最適推定値算出部１６８において、時刻ｎまでの情報による時刻ｎ＋１での所望信号の最適推定値ベクトルｘｐ１（ｎ＋１｜ｎ）を計算する。この計算は、ステップＳ１０００で設定した状態遷移行列Φｐ１、および、前回のステップＳ１４００で計算した所望信号の最適推定値ベクトルｘｐ１（ｎ｜ｎ）を用いて、次の式（８）により行う。なお、このステップＳ１３００は、図１２の反復過程の手順３に相当する。
【００７６】
【数８】
0005447867.tif 000009
【００７７】
次に、同じく最適推定値算出部１６８で、（ｎ＋１）→（ｎ＋１）の状態量（所望信号）の最適推定値（ベクトル）を計算する（Ｓ１４００）。具体的には、最適推定値算出部１６８において、時刻ｎ＋１までの情報による時刻ｎ＋１での所望信号の最適推定値ベクトルｘｐ１（ｎ＋１｜ｎ＋１）を計算する。この計算は、ステップＳ１３００で計算した所望信号の最適推定値ベクトルｘｐ１（ｎ＋１｜ｎ）、ステップＳ１２００で計算した重み係数行列Ｋｐ１（ｎ＋１）、ステップＳ１０００で設定した観測遷移行列Ｍｐ１、および時刻ｎ＋１における観測信号ベクトルｙｐ１（ｎ＋１）を用いて、次の式（９）により行う。なお、このステップＳ１４００は、図１２の反復過程の手順４に相当する。
【００７８】
【数９】
0005447867.tif 000010
【００７９】
次に、処理を終了するか否かを判断する（Ｓ１５００）。この判断は、例えば、時刻ｎが所定のサンプル数Ｎに達したか否かを判定することによって行う。この判断の結果として、時刻ｎが所定のサンプル数Ｎに達していない場合は（Ｓ１５００：ＮＯ）、ステップＳ１６００に進み、時刻ｎが所定のサンプル数Ｎに達した場合は（Ｓ１５００：ＹＥＳ）、ステップＳ１８００に進む。なお、判断の基準は、上記の例に限定されない。例えば、リアルタイムで処理を行う場合は、時刻ｎが所定のサンプル数Ｎに達していなくても、サンプルがなくなった時点で処理を終了するようにしてもよい。
【００８０】
ステップＳ１６００では、相関演算部１６４で、（ｎ＋１）→（ｎ＋１）の誤差共分散行列、つまり、（ｎ＋１）→（ｎ＋１）の推定誤差の相関値（行列）を計算する。具体的には、相関演算部１６４において、時刻ｎ＋１までの情報により時刻ｎ＋１での所望信号ベクトルを推定した場合の誤差（所望信号の推定誤差ベクトル）の相関行列Ｐｐ１（ｎ＋１｜ｎ＋１）を計算する。この計算は、ステップＳ１２００で計算した重み係数行列Ｋｐ１（ｎ＋１）、ステップＳ１０００で設定した観測遷移行列Ｍｐ１、およびステップＳ１１００で計算した所望信号の推定誤差ベクトルの相関行列Ｐｐ１（ｎ＋１｜ｎ）を用いて、次の式（１０）により行う。なお、このステップＳ１６００は、図１２の反復過程の手順５に相当する。
【００８１】
【数１０】
0005447867.tif 000011
【００８２】
次に、ステップＳ１７００では、時刻ｎのカウンタを１だけインクリメントして（ｎ＝ｎ＋１）、ステップＳ１１００に戻る。
【００８３】
一方、ステップＳ１８００では、本アルゴリズムの計算結果を出力値として一時保存する。具体的には、ステップＳ１４００で計算した所望信号の最適推定値ベクトルｘｐ１（ｎ＋１｜ｎ＋１）を、本アルゴリズムの出力値として画像復元処理部１６０内または記憶部１５０内に一時保存する。
【００８４】
図１４は、発明手法を視覚的にまとめた説明図である。このように、発明手法では、新しい状態空間モデル（状態方程式と観測方程式）を構成しているため、１段階処理で画像の復元処理が可能となる。このことは、本発明の大きな特徴の１つである。
【００８５】
ここで、カルマンフィルタを用いた従来の画像復元手法との対比において、発明手法の特徴と効果をまとめて説明しておく。
【００８６】
カルマンフィルタを用いた従来の画像復元手法では、上記のように、２段階処理（ステップ１で、ＡＲ次数を決定し、ＡＲ係数を推定した後、ステップ２で、この推定したＡＲ係数を用いて状態空間モデル（状態方程式と観測方程式）を構成してカルマンフィルタを実行する）によって、画像の復元を実現している。したがって、ステップ１でのＡＲ次数の決定およびＡＲ係数の推定の精度によってステップ２でのカルマンフィルタによる画像復元能力が大きく左右されることが、容易に予想される。これに対して、発明手法では、ＡＲシステムのコンセプトを必要としない新しい状態空間モデル（状態方程式と観測方程式）を構成し、これを用いた１段階処理の新しい予測法によって、高性能な画像復元を実現している。また、発明手法では、処理ステップを１段階減らすことができるため、演算量を削減することができ、ひいては、回路規模の縮小やメモリ容量の削減などを図ることができる。
【００８７】
また、カルマンフィルタを用いた従来の画像復元手法では、ステップ１でＡＲ係数の推定を行う際に必要となるＡＲ係数の次数の決定が大きな問題点となる。一般にＡＲ係数の次数は状態量に依存しているため、状態量が既知でなければ、理論的にＡＲ係数の次数を正確に決定することは困難である。このことは、状態量が既知でなければならないことを意味しているため、リアルタイム処理は困難となる。また、それにより、正確でないＡＲ係数の次数を用いることになるため、正確なＡＲ係数の推定は困難である。したがって、このことは、カルマンフィルタを用いた従来の画像復元手法の画像復元能力を低下させる大きな要因となっている。また、何らかの技術でＡＲ次数とＡＲ係数をリアルタイムに推定できるようになったとしても、処理ステップが１段階増えることより演算量の増加は避けられない。これに対して、発明手法では、ＡＲシステムのコンセプトを必要としないため、このような問題は生じない。
【００８８】
また、カルマンフィルタを用いた従来の画像復元手法では、状態量をＡＲシステムを用いて表現することによってモデル化している。このことは、カルマンフィルタを用いた従来の画像復元手法はＡＲシステムでモデル化可能な状態量にしか適用できないことを意味している。すなわち、カルマンフィルタを用いた従来の画像復元手法は、ＡＲシステムでモデル化が困難な状態量には適用することができない。これに対して、発明手法は、ＡＲシステムのコンセプトを必要としないため、適用対象に関するこのような制約は存在しない。
【００８９】
また、カルマンフィルタを用いた従来の画像復元手法では、状態方程式の駆動源が白色信号で、かつ、状態量と観測雑音が無相関であることを仮定した上で、カルマンフィルタ理論を適用している。これに対して、発明手法では、駆動源が有色信号（クリアな画像）となるものの、状態方程式と観測方程式の特殊な構成によって発明手法のアルゴリズムを実行することができる。このことは、発明手法は、カルマンフィルタ理論の一般的な適用条件を考慮しなくても実行可能であることを意味している。すなわち、発明手法は、カルマンフィルタ理論よりも適用範囲が広いといえる。
【００９０】
したがって、ＡＲシステムのコンセプトを必要としない発明手法は、画像の撮り直しが許されない場面、例えば、医療分野における瞬間的な心臓や肺などの画像復元、汚れなどで劣化した古文書などの復元、および文字・物体認識においても貢献可能な技術である。
【００９１】
また、本発明者は、本発明の効果（発明手法の有効性）を実証するためにシミュレーションを行った。具体的には、発明手法１の画像復元能力を評価するために、（１）視覚評価、（２）客観評価、および（３）主観評価を行った。視覚評価は、原画像と復元画像とを視覚により比較した評価である。客観評価は、数値による評価である。主観評価は、聞き取り調査である。以下、これらを順に説明する。
【００９２】
まず、シミュレーション条件について説明しておく。
【００９３】
図１５は、シミュレーション条件を説明するための図である。
【００９４】
本シミュレーションでは、原画像として、図１５に示す２つの画像、つまり、（ａ）「カメラマン（Cameraman）」と（ｂ）「レナ（Lenna）」を用いた。また、原画像に加えるぼけに対応する点拡がり関数（ＰＳＦ）のモデルとして、図１５に示す２次元ガウス関数を用い、雑音として、加法性白色ガウス雑音を用いた。また、信号対雑音比ＳＮＲを３０ｄＢとした。また、比較のため、シミュレーションは、同一の上記シミュレーション条件の下で、従来手法１（ウィナーフィルタを用いた画像復元手法）と、従来手法２（射影フィルタを用いた画像復元手法）と、従来手法３（カルマンフィルタを用いた画像復元手法）と、発明手法（発明手法１）に対して行った。例えば、画像「カメラマン」は、図１６に示すように、空の部分は画像の分散が変化していないため定常性が強く、人物の頭の部分は画像の分散が変化しているため非定常性が強いといえる。
【００９５】
（１）視覚評価
図１７は、原画像「カメラマン」に対するシミュレーション結果（視覚評価）を示す図であり、図１８は、図１７中の破線の円の周辺領域を拡大した図である。
【００９６】
図１８によく示すように、従来手法１は、ぼけがほとんど取り除かれていない。また、従来手法２は、従来手法１よりもぼけが取り除かれているが、原画像と比べるとまだぼけが残っている。一方、従来手法３は、従来手法１、２に比べて、ぼけは取り除かれてくっきりしているように見えるが、原画像とは異なった復元画像となっている。特に、従来手法３は、原画像と比べて全体的に暗くなっており、空の部分においては、劣化画像よりも劣化してしまっている。
【００９７】
これに対して、発明手法は、図１８によく示すように、従来手法１、２、３に比べて原画像を忠実に復元していることが確認できる。すなわち、図１８の拡大画像においても、発明手法の有効性を確認することができる。
【００９８】
図１９は、原画像「カメラマン」に対するシミュレーション結果（視覚評価）を示す図であり、図２０は、図１９中の破線の円の周辺領域を拡大した図である。
【００９９】
図２０によく示すように、カメラの三脚部分に注目すると、従来手法１は、全くぼけが取り除かれていない。また、従来手法２は、従来手法１よりはぼけが取り除かれているが、原画像までは復元されていない。一方、従来手法３は、従来手法１、２に比べるとぼけは取り除かれているように見えるが、画像全体が荒れてしまっており、原画像とかけ離れた復元画像となっていることが確認できる。
【０１００】
これに対して、発明手法は、図２０によく示すように、従来手法１、２、３に比べて原画像をより忠実に復元していることが確認できる。すなわち、図２０の拡大画像においても、発明手法の有効性を確認することができる。
【０１０１】
図２１は、原画像「レナ」に対するシミュレーション結果（視覚評価）を示す図であり、図２２は、図２１中の破線の円の周辺領域を拡大した図である。
【０１０２】
図２２によく示すように、従来手法１は、ほとんどぼけが取り除かれておらず、復元画像が原画像と比べて全体的に白くなっている（輝度が高くなっている）。また、従来手法２は、従来手法１に比べるとぼけは取り除かれているが、原画像まではぼけが取り除かれていない。一方、従来手法３は、従来手法１、２に比べるとぼけが取り除かれておりはっきりとしているように見えるが、劣化画像よりも劣化してしまっており、原画像とはかけ離れた画像になっている。特に肌の部分に注目すると、この点が顕著に表れている。
【０１０３】
これに対して、発明手法は、図２２によく示すように、従来手法１、２、３に比べて原画像を忠実に復元していることが確認できる。すなわち、図２２の拡大画像においても、発明手法の有効性を確認することができる。
【０１０４】
（２）客観評価（数値による評価）
図２３は、原画像に対するシミュレーション結果（客観評価）を示す図である。
【０１０５】
ここでは、従来手法と発明手法の画像復元能力を数値により評価するため、図２３にも示す次の式（１１）で表されるＳＮＲｏｕｔ［ｄＢ］を用いて画像復元能力を評価した。なお、ＳＮＲは、信号に対する雑音の比であり、数値が大きいほど劣化（ぼけや雑音など）が少なく画像として良好であるといえる。
【０１０６】
【数１１】
0005447867.tif 000012
【０１０７】
この場合、図２３に示すように、「カメラマン」、「レナ」どちらの画像についても、発明手法は、従来手法１、２、３よりもＳＮＲｏｕｔの数値が大きいことが確認できる。これにより、発明手法は、客観評価の点でも、従来手法１、２、３に比べて画像復元能力が高いことがわかる。
【０１０８】
（３）主観評価（聞き取り調査）
図２４は、原画像に対するシミュレーション結果（主観評価）を示す図である。
【０１０９】
ここでは、従来手法と発明手法の画像復元能力を評価するために、聞き取り調査による主観評価を行った。画像復元の性能評価は、ＡＣＲ（絶対範疇評価）に基づく５段階ＭＯＳ（平均オピニオン値）を用いた聞き取り調査により行った。ＭＯＳ評価基準は、図２４に示す通りである。５０人の聴取者が画像復元により得られた画像（図１７〜図２２参照）を評価した。各々の聴取者は、評価値「１」から「５」を決定する。評価値「５」が最良である。
【０１１０】
図２４に示すように、ＭＯＳ評価は、「カメラマン」、「レナ」どちらの画像についても、発明手法は、従来手法１、２、３よりも高い評価を得ていることが確認できる。これにより、発明手法は、主観評価の点でも、従来手法１、２、３に比べて画像復元能力が高いことがわかる。
【０１１１】
以上のシミュレーション結果により、本発明の画像復元方法（発明手法１）は、従来手法に比べて高い画像復元能力を発揮していることがわかる。特に、図１７〜図２２の視覚評価から明らかなように、本発明の画像復元方法は、非定常性が強いエッジ部分などにおいて、従来手法に比べて著しく復元精度が高いことがわかる。
【０１１２】
次いで、上記構成を有する画像復元装置１００の動作について、図２５に示すフローチャートを用いて説明する。なお、図２５に示すフローチャートは、記憶部１５０（主記憶装置１５２または補助記憶装置１５４）に制御プログラムとして格納されており、図示しないＣＰＵによって実行される。
【０１１３】
まず、ステップＳ２０００では、画像入力装置１１０から復元処理対象の画像データ（劣化画像）を取り込み、記憶部１５０（主記憶装置１５２または補助記憶装置１５４）の所定の記憶領域に格納する。
【０１１４】
そして、ステップＳ２０５０では、ステップＳ２０００で取り込んだ画像データを表示メモリ１５６に書き込んで、ディスプレイ１８２に表示する。
【０１１５】
そして、ステップＳ２１００では、操作部１３０の領域指定部１３４で、領域指定処理を行う。具体的には、使用者の入力操作により、ステップＳ２０５０でディスプレイ１８２に表示した画像のうち画像復元処理の対象となる領域（画像の特定の範囲）が指定されると、その指定された領域のデータを生成する。領域の指定は、例えば、画面上、ポインタなどで行われる。なお、使用者により領域の指定が行われなかった場合は、表示された画像全体を指定したものとして取り扱われる。
【０１１６】
そして、ステップＳ２１５０では、拡大処理を行うか否かを判断する。この判断は、例えば、ステップＳ２１００で生成した指定領域データに基づいて行われる。具体的には、指定された領域が表示された画像全体よりも小さい場合には、拡大処理を行うと判断する。この判断の結果として、拡大処理を行う場合は（Ｓ２１５０：ＹＥＳ）、ステップＳ２２００に進み、拡大処理を行わない場合は（Ｓ２１５０：ＮＯ）、直ちにステップＳ２２５０に進む。
【０１１７】
ステップＳ２２００では、ステップＳ２１００で指定した領域の拡大処理を行う。具体的には、例えば、指定された領域が表示された画像全体に対応する大きさとなるように拡大処理を行う。拡大処理の結果は、記憶部１５０（主記憶装置１５２または補助記憶装置１５４）の作業用メモリに書き込まれる。なお、拡大処理が行われなかった場合も、その画像データは作業用メモリに書き込まれる。
【０１１８】
そして、ステップＳ２２５０では、時刻ｎの処理対象ブロック（例えば、３×３の大きさ）を選択する。
【０１１９】
そして、ステップＳ２３００では、ステップＳ２２５０で選択した処理対象ブロックに対して、画像の急激な変化領域が検出されたか否かを判断する。画像の急激な変化領域は、例えば、画像のエッジ部分などに相当する。具体的には、例えば、作業用メモリに書き込まれた画像データ（拡大処理の有無を問わず）に対して、ステップＳ２２５０で選択した処理対象ブロック（３×３の大きさ）上で画素データを順次スキャンすることにより、急激な画素データの変化点があるか否かを検出する。この判断の結果として、画像の急激な変化領域が検出された場合は（Ｓ２３００：ＹＥＳ）、ステップＳ２３５０に進み、画像の急激な変化領域が検出されなかった場合は（Ｓ２３００：ＮＯ）、ステップＳ２４００に進む。
【０１２０】
具体的な検出方法の一例は、次の通りである。ある一定量の過去の観測量の平均値ｙ’（ｎ−１）と現在の観測量ｙ（ｎ）との差の絶対値、つまり、|ｙ’（ｎ−１）−ｙ（ｎ）|の値をしきい値αと比較する。そして、その値がしきい値α以上の場合、つまり、|ｙ’（ｎ−１）−ｙ（ｎ）|≧αの場合は、急激な画素データの変化点があると判断し、その値がしきい値α未満の場合、つまり、|ｙ’（ｎ−１）−ｙ（ｎ）|＜αの場合は、急激な画素データの変化点はないと判断する。例えば、図２６において、時刻ｎ−１の処理対象ブロック３０１における所定部分３０２の観測量の平均値ｙ’（ｎ−１）と、時刻ｎの処理対象ブロック３０３における注目画素の部分３０４の観測量ｙ（ｎ）とを比較する。
【０１２１】
ステップＳ２３５０では、ステップＳ２２５０で選択した処理対象ブロックに対して、発明手法による画像復元処理を行う。発明手法による画像復元処理は、上記のように、急激な画素データの変化点がある場合、つまり、エッジ部分などを含む場合であっても、高精度な復元処理が可能である。なお、発明手法による画像復元処理の詳細な手順の一例については、図１３のフローチャートを用いて説明した通りである。
【０１２２】
一方、ステップＳ２４００では、ステップＳ２２５０で選択した処理対象ブロックに対して、発明手法以外の他の画像復元手法による画像復元処理を行う。他の画像復元手法としては、ウィナーフィルタを用いた画像復元手法や射影フィルタを用いた画像復元手法、カルマンフィルタを用いた画像復元手法などを含む任意の画像復元手法を用いることができる。急激な画素データの変化点がない場合、つまり、エッジ部分などを含まない場合には、他の画像復元手法でも精度良く復元処理を行うことができるからである。
【０１２３】
そして、ステップＳ２４５０では、ステップＳ２３５０での復元処理結果またはステップＳ２４００での復元処理結果を、順次、記憶部１５０（主記憶装置１５２または補助記憶装置１５４）の作業用メモリに保存する。
【０１２４】
そして、ステップＳ２５００では、時刻ｎのカウンタの値を１だけインクリメントする。
【０１２５】
そして、ステップＳ２５５０では、１頁分の画像復元処理が終了したか否かを判断する。この判断は、時刻ｎのカウンタの値に基づいて行われる。この判断の結果として、１頁分の画像復元処理が終了した場合は（Ｓ２５５０：ＹＥＳ）、ステップＳ２６００に進み、１頁分の画像復元処理が終了していない場合は（Ｓ２５５０：ＮＯ）、ステップＳ２２５０に戻る。
【０１２６】
ステップＳ２６００では、１頁分の画像復元処理が終了したため、表示メモリ１５６の更新処理を行う。すなわち、ここでは、１頁分の画像復元処理が終了した時点で表示メモリ１５６を更新するように構成されている。
【０１２７】
そして、ステップＳ２６５０では、ステップＳ２０００で取り込んだ画像データの全頁分に対して画像復元処理が終了したか否かを判断する。この判断の結果として、全頁分の画像復元処理が終了していない場合は（Ｓ２６５０：ＮＯ）、ステップＳ２２５０に戻り、全頁分の画像復元処理が終了した場合は（Ｓ２６５０：ＹＥＳ）、上記一連の処理を終了する。
【０１２８】
次に、ぼけ情報である点拡がり関数（ＰＳＦ）に対応する観測遷移行列Ｍｐ１の設定と実装の具体例について、図２７および図２８を用いて説明する。
【０１２９】
図２７は、復元モードの実装の一例を説明するための図である。
【０１３０】
上記のように、観測遷移行列Ｍｐ１は、式（４）で定義された９×９次の行列であって、画像の劣化モデルにおけるぼけの点拡がり関数（ＰＳＦ）に対応しており、観測遷移行列Ｍｐ１を構成する各要素ｈｒ，ｓ（ｒ，ｓ＝−１，０，１）は、既知であって、予めデータ化して適切に設定されている。すなわち、本実施の形態では、画像によってぼけ具合が変わるため、いわばぼけ情報である観測遷移行列Ｍｐ１の各要素ｈｒ，ｓの値の組を復元モードとしていくつか予め設定しておき、使用者は、操作部１３０の復元モード指定部１３６を介して、その復元モードを任意に指定できるように構成されている。例えば、図２７の例では、デジタルカメラなどによく設定されている夜景モードやスポーツモードなどの撮影モード毎に、予め、その撮影モードに合った適切な観測遷移行列Ｍｐ１がそれぞれ設定されている。したがって、使用者は、撮影モードを切り替えることによって、最適な復元モードを使用することができる。なお、復元モード（観測遷移行列Ｍｐ１の各要素ｈｒ，ｓの値）自体を自動または手動で再調整することも可能である。
【０１３１】
図２８は、復元モードの実装の他の例を説明するための図である。
【０１３２】
図２８の例では、画像によってぼけ具合が変わるため、まず、図２８（Ａ）に示すように、ディスプレイ１８２上で、処理対象の劣化画像を選択した後、図２８（Ｂ）に示すように、観測遷移行列Ｍｐ１（各要素ｈｒ，ｓの値）を変えた復元画像をいくつかディスプレイ１８２に表示する。このとき、使用者が、例えば、操作部１３０（タッチパネルなど）を介して、ディスプレイ１８２に表示された複数の復元画像の中から、最適な復元画像を選択すると、復元画像が決定される。なお、使用者の選択の履歴をデータとして蓄積することによって、使用者の好みに合った復元モード（観測遷移行列Ｍｐ１の各要素ｈｒ，ｓの値）を自動で選択することも可能となる。
【０１３３】
このように、本実施の形態によれば、ＡＲシステムのコンセプトを必要としない新しい状態空間モデル（状態方程式と観測方程式）を構成して１段階処理の新しい予測法によって画像復元を行うため、シンプルで実用的かつ高性能な画像復元を実現することができる。すなわち、ＡＲ次数の決定およびＡＲ係数を推定するステップを必要としないシンプルな構成で、かつ、非定常性が強い自然画像に対しても有効に復元処理を行うことができる実用性を有し、しかも、従来手法に比べて画像復元能力を向上することができる。
【０１３４】
また、本実施の形態では、上記のように、ＡＲシステムのコンセプトを必要としないため、画像の撮り直しが許されない場面にも広く適用可能である。なお、もちろん、ここでも、上記のように、画像の種類は限定されず、静止画像でも動画像でもよい。
【０１３５】
例えば、古文書の傷や汚れなどをぼけや雑音とみなし、発明手法を用いてこれらを取り除き、文字・物体認識技術を用いることによって、古文書の解読、復元、およびデジタルデータベース化に効果を発揮することが期待される。
【０１３６】
また、防犯分野においても、近年、防犯カメラ（監視カメラ）などの設置が普及しつつあり、犯罪者の多くは防犯カメラなどに撮影されていることが多いため、この画像が犯罪者の発見の手がかりとなることがある。しかし、これらの画像は一般に画質が荒く、動いている対象物の画像はぼけが生じていることが一般的である。また、これらの画像を拡大した場合には、通常、画像はより一層劣化する。そこで、拡大画像を含めてこれらの劣化画像に発明手法を適用することによって、劣化画像のぼけと雑音を取り除き、クリアな画像を提供することで、犯罪者の早期発見につながることが期待される。なお、このように、拡大画像を含めて劣化画像に発明手法を適用することによって、劣化画像のぼけと雑音を取り除き、クリアな画像を提供することは、防犯分野に限らず、監視カメラなどで撮影された画像に基づいて、事故の判定や器具の故障診断などを行う場合にも適用可能である。
【０１３７】
また、近年、カメラ機能付き携帯電話やデジタルカメラなどが急速な勢いで普及し、１人１台の時代となっている。しかし、これらの製品に搭載されている画像復元技術は、顔認識やフィルタなどのぼけや雑音を予防する技術であり、撮り直しが許されない場面を想定していない。したがって、これらの製品に発明手法を適用することによって、撮り直しが許されない瞬間的な画像復元が可能になる。
【０１３８】
また、医療の分野において、健康を調査する最も効果的な手段の１つとして、内視鏡やレントゲン、エコー、ＣＴ、ＭＲＩなどの画像診断がある。例えば、肺や心臓などの画像診断では、肺や心臓などの動作によってぼけと雑音を含んだ情報から健康診断を行っている。これらの画像診断に発明手法を適用してぼけと雑音を取り除くことによって、専門医にクリアな画像を迅速に提供することが可能となり、病気の早期発見の一助になることが期待される。
【０１３９】
また、近年カーナビゲーションシステムが普及していることから、フロント・リアカメラを搭載した自動車も増えている。このようなリアルタイムで画像復元処理が必要な自動車のフロント・リアカメラなどの場合においても発明手法は有効である。
【０１４０】
２００８年８月８日出願の特願２００８−２０６３１６の日本出願に含まれる明細書、図面および要約書の開示内容は、すべて本願に援用される。
【産業上の利用可能性】
【０１４１】
本発明に係る画像復元装置および画像復元方法は、画像復元能力を向上することができるシンプルで実用的な画像復元装置および画像復元方法として有用である。
【符号の説明】
【０１４２】
１００  画像復元装置
１１０  画像入力装置
１１２  カメラ
１１４  スキャナ
１１６、１８６  記録メディア
１１８、１８８  モデム
１２０  入力インタフェース部
１３０  操作部
１３２  パラメータ設定部
１３４  領域指定部
１３６  復元モード指定部
１４０  内部インタフェース部
１５０  記憶部
１５２  主記憶装置
１５４  補助記憶装置
１５６  表示メモリ
１６０  画像復元処理部
１６０ａ  第１復元処理部
１６０ｂ  第２復元処理部
１６２  初期設定部
１６４  相関演算部
１６６  重み係数算出部
１６８  最適推定値算出部
１７０  出力インタフェース部
１８０  画像出力装置
１８２  ディスプレイ
１８４  プリンタ
【図１】
0005447867.tif 000013
【図２】
0005447867.tif 000014
【図４】
0005447867.tif 000015
【図５】
0005447867.tif 000016
【図６】
0005447867.tif 000017
【図７】
0005447867.tif 000018
【図８】
0005447867.tif 000019
【図９】
0005447867.tif 000020
【図１０】
0005447867.tif 000021
【図１１】
0005447867.tif 000022
【図１２】
0005447867.tif 000023
【図１３】
0005447867.tif 000024
【図１４】
0005447867.tif 000025
【図２３】
0005447867.tif 000026
【図２４】
0005447867.tif 000027
【図２５】
0005447867.tif 000028
【図２６】
0005447867.tif 000029
【図２７】
0005447867.tif 000030
【図２８】
0005447867.tif 000031
【図３】
0005447867.tif 000032
【図１５】
0005447867.tif 000033
【図１６】
0005447867.tif 000034
【図１７】
0005447867.tif 000035
【図１８】
0005447867.tif 000036
【図１９】
0005447867.tif 000037
【図２０】
0005447867.tif 000038
【図２１】
0005447867.tif 000039
【図２２】
0005447867.tif 000040