説明

画像復元装置および画像復元方法

【課題】画像復元性能を向上できるシンプルで実用的な画像復元装置を提供する。
【解決手段】画像復元装置は、時刻nのみの劣化画像情報に対して、時刻nまたは時刻n+1までの情報により原画像情報を含む時刻n+1でのシステムの状態量を推定した場合の推定誤差の相関値を算出する相関演算部と、時刻nのみの劣化画像情報に対して、前記相関演算部によって算出された相関値を用いて、時刻n+1までの情報による時刻n+1での前記状態量の最適推定値と、時刻nまでの情報による時刻n+1での前記状態量の最適推定値と、前記劣化画像情報を含む観測量の推定誤差と、の関係を規定するための重み係数を算出する重み係数算出部と、時刻nのみの劣化画像情報に対して、前記重み係数算出部によって算出された重み係数を用いて、時刻nまたは時刻n+1までの情報による時刻n+1での前記状態量の最適推定値を算出する最適推定値算出部とを有する。

【発明の詳細な説明】
【技術分野】
【0001】
本発明は、画像復元装置および画像復元方法に関する。
【背景技術】
【0002】
近年、画像工学の分野において、劣化した画像から原画像を復元する技術の研究開発が多くなされている。すなわち、原画像(所望情報、クリアな画像)に不必要な情報(ぼけと雑音)が混在した劣化画像(受信情報)から不必要な情報(ぼけと雑音)を取り除き、原画像(所望情報)のみを抽出することは、画像復元の分野において必要不可欠な技術であり、近年盛んに研究開発が行われている。例えば、デジタルカメラ(デジタルスチルカメラおよびデジタルビデオカメラの総称)や携帯電話などで撮影された画像は、手ぶれや焦点ずれなどに起因する「ぼけ」と、暗電流や熱雑音などに起因するガウス性またはインパルス性の「雑音」との影響から、実物と比較して画質が劣化することは避けられない。「画像復元」とは、このように劣化した画像からできるだけ原画像に近い画像を復元することである。
【0003】
現在一般に市場で普及している画像復元技術の多くは、例えば、手ぶれ補正や顔認識、色補正、様々なフィルタなどを用いて、予めぼけと雑音の影響を軽減させる予防的な技術が主流となっている。この結果、近年、特にデジタルカメラの分野では、デジタルカメラの多機能化・高性能化により鮮明な画像が気軽に手に入れられるようになってきている。
【0004】
しかし、これらの予防的な技術は、何度も画像を撮り直せる状況下では問題ないが、古文書などの既に劣化している画像や、スポーツや医療などの分野における瞬間的に変化する画像など、撮り直しが許されない画像に対する復元は、いまだ困難な問題である。ここで、スポーツと医療の分野における瞬間的に変化する画像は、例えば、選手の瞬間的な動作や、肺や心臓などの臓器の瞬間的な様子などである。したがって、今日、特に、撮り直しが許されない状況下における画像復元の重要性が高まってきている。
【0005】
撮り直しが許されない状況下における従来の画像復元技術として、ウィナーフィルタ(Wiener Filter)を用いた画像復元手法が広く知られている(非特許文献1、非特許文献2)。この手法は、フィルタを通して得られる復元画像と原画像との平均2乗誤差を最小にするフィルタであって、最小2乗フィルタとも呼ばれている。この手法は、周波数領域で処理を行うため、確率過程の定常性と半無限長の画像のサイズとを前提にした画像復元手法である。
【0006】
また、他の画像復元技術として、射影フィルタを用いた画像復元手法も知られている(非特許文献3、非特許文献4)。射影フィルタは、原画像と復元画像との近さを直接原画像の空間で評価しており、復元画像の雑音の画像成分が個々の原画像の最良近似、つまり、個々の原画像の正射影となるものの中で、復元画像の雑音成分の2乗平均を最小にするものである。また、この性質から、射影フィルタは、出現頻度に無関係に、最良近似画像を復元する手法である。
【0007】
また、さらに他の画像復元技術として、カルマンフィルタ(Kalman Filter)を用いた画像復元手法も知られている(非特許文献5、非特許文献6)。この手法は、まずステップ1で、劣化画像を用いて原画像のAR(Auto Regressive:自己回帰)システムのパラメータ(以下「AR係数」という)を推定した後、ステップ2で、ステップ1で推定したAR係数を用いて状態空間モデル(状態方程式と観測方程式)を構成し、これにカルマンフィルタ理論(カルマンフィルタアルゴリズム)を適用することによって、高性能な画像復元を実現している。さらには、このような2段階処理を行うカルマンフィルタを2次元カルマンフィルタに拡張した画像復元手法も知られている(非特許文献7、非特許文献8)。この手法は、2次元に拡張した状態空間モデルを用いて、劣化画像を復元するカルマンフィルタを実現している。なお、非特許文献7および非特許文献8に示す2次元カルマンフィルタと区別するため、適宜、非特許文献5および非特許文献6に示すカルマンフィルタを「1次元カルマンフィルタ」と称する。
【先行技術文献】
【非特許文献】
【0008】
【非特許文献1】西宮亮平、外1名、「複数のウィーナーフィルタによる画像復元」(Image Restoration by Using Multiple Wiener Filters)、電子情報通信学会論文誌、A、Vol.J83−A、No.7、pp.892−902、2000年7月
【非特許文献2】山根延元、外3名、「画像の局所定常ガウス混合モデルに基づく適応的ウィーナーフィルタを用いた最適雑音除去」(An Optimal Noise Removal Using Adaptive Wiener Filter Based on Locally Stationary Gaussian Mixture Distribution Model for Images)、電子情報通信学会論文誌、A、Vol.J85−A、No.9、pp.993−1004、2002年9月
【非特許文献3】小川 英光 、外1名、「部分射影フィルタの諸性質」(Properties of Partial Projection Filter)、電子情報通信学会論文誌、A、Vol.J71−A、No.2、pp.527−534、1988年2月
【非特許文献4】小出 裕司、外2名、「信号・画像推定のための射影フィルタ族の統一理論」(A Unified Theory of the Family of Projection Filters for Signal and Image Estimation)、電子情報通信学会論文誌、D−II、Vol.J77−D−II、No.7、pp.1293−1301、1994年7月
【非特許文献5】城崇、外2名、「カルマンフィルタによる高精度画像復元のための画像のモデリングとパラメータ推定」(Image Modeling and Parameter Identification for Image Restoration Using a Kalman Filter)、電子情報通信学会論文誌、D−II、Vol.J80−D−II、No.11、pp.2912−2919、1997年11月
【非特許文献6】松村淳、外2名、「効果的な雑音除去のための適応的な画像のモデル化によるカルマンフィルタ」(A Kalman Filter Using Adaptive Image Modeling for Noise Reduction)、電子情報通信学会論文誌、D−II、Vol.J86−D−II、No.2、pp.212−222、2003年2月
【非特許文献7】Mahmood R Azimi-Sadjadi and Ping Wong, "Two-Dimensional Block Kalman Filtering for Image Restoration," IEEE Trans. On Signal Process, vol.ASSP-35, no.12, Dec. 1987
【非特許文献8】Stuart Citrin and Mahmood R. Azimi-Sadjadi, "A Full-Plane Block Kalman Filter for Image Restoration," IEEE Trans. on Image Processing, vol.1, no.4, Oct. 1992
【発明の概要】
【発明が解決しようとする課題】
【0009】
しかしながら、ウィナーフィルタを用いた画像復元手法は、画像の劣化状態によらず復元が可能であるという利点を有する反面、非定常性(画像の分散が変化していること)が強い自然画像(加工されていない撮影されたままの画像)に対する復元精度が低いという欠点がある。
【0010】
すなわち、ウィナーフィルタを用いた画像復元手法は、上記のように、周波数領域で処理を行うため、確率過程の定常性と半無限長の画像のサイズとを前提にしているが、現実には、実環境においてこの前提条件が成り立つことは困難であるため、適切に復元されない場合が存在する(自然画像の非定常性の問題)。また、この手法は、最小2乗誤差を評価量とした一括処理であるため、復元画像にぼけが残るという問題が生じる可能性もある(評価量の問題)。
【0011】
具体的には、ウィナーフィルタを用いた画像復元手法では、画像の出現頻度が低くなると、復元精度が低下する、つまり、画像においては、出現頻度が低いエッジ部分の復元に影響するという問題がある。換言すれば、例えば、カメラで撮影された画像は、エッジ部分が多く非定常性が強いため、特にエッジ部分において最適に復元できないという問題がある。ただし、階調や色などが変化していない、つまり、画像の分散が変化していない、単調な背景部分などは、定常性が強いため、通常、良好な復元精度を得ることができる。
【0012】
また、カルマンフィルタを用いた画像復元手法は、ウィナーフィルタを用いた画像復元手法の問題点(自然画像の非定常性の問題と評価量の問題)を解決するための手法であって、ウィナーフィルタを用いた画像復元手法の問題点を解決可能であるという利点を有する反面、処理対象画像(観測画像・劣化画像)にぼけに基づく劣化があると、復元精度が低いという問題がある。
【0013】
すなわち、カルマンフィルタを用いた画像復元手法は、ステップ1で、注目画素とその周囲の画素との相関を考慮すべく、前処理として注目画素とその周囲の画素に対して原画像のAR係数を推定した後、ステップ2で、ステップ1で推定したAR係数を用いて状態空間モデル(状態方程式と観測方程式)を構成し、具体的には、ステップ1で推定したAR係数と原画像から状態方程式を構成し、原画像とぼけ関数と雑音から観測方程式を構成することにより、カルマンフィルタ理論(カルマンフィルタアルゴリズム)を用いて画像を復元している。したがって、カルマンフィルタを用いた画像復元手法は、定常性を前提としない時間領域のみの処理であり、かつ、推定誤差の分散を評価量とする逐次的処理であるため、ウィナーフィルタを用いた画像復元手法の問題点を解決することができる。なお、このような2段階処理を必要とする点は、1次元カルマンフィルタを用いた手法も2次元カルマンフィルタを用いた手法も同様である。
【0014】
しかしながら、一方で、カルマンフィルタを用いた手法は、1次元カルマンフィルタであれ2次元カルマンフィルタであれ、ステップ1で推定したAR係数を用いてステップ2のカルマンフィルタアルゴリズムを実行するため、ステップ1でのAR係数の推定精度に劣化画像の復元精度が大きく依存してしまうという問題を有している(ARシステムの問題)。例えば、デジタルカメラの場合、処理対象画像にぼけ(焦点ずれなど)に基づく劣化が含まれると、ステップ1でのAR次数の決定およびAR係数の推定が困難となるため、ステップ2でのカルマンフィルタによる画像復元の精度に影響が出てしまう。
【0015】
この点、一般には、AR係数の正確な推定は困難である。AR係数の正確な推定は、例えば、画像復元であれば、クリアな画像(原画像)に依存しているからである。このことは、原画像が既知でなければならないことを意味しているため、リアルタイム処理は困難となる。また、仮に何らかの手法でリアルタイムにAR係数を正確に推定することが可能になったとしても、処理が増加するため演算量の問題は避けられない。また、そもそも、AR係数の推定はAR係数の次数を決定した後に行われるが、AR係数の次数を決定するのは非常に困難であり、この点からもAR係数の正確な推定は困難であるといえる。
【0016】
したがって、今日、ウィナーフィルタを用いた画像復元手法の問題点およびカルマンフィルタを用いた画像復元手法の問題点の両方を解決しうる、シンプルで実用的な高性能の画像復元手法、つまり、シンプルな構成を有し、かつ、実環境でも使用可能であり、しかも劣化画像の復元性能(つまり、画像復元能力)が高い画像復元手法が求められている。特に、劣化画像の復元性能(画像復元能力)としては、クリアな画像にぼけと雑音が含まれた劣化画像に対して、ぼけと雑音を同時に取り除くことができる高い性能(能力)が求められている。
【0017】
本発明の目的は、ぼけと雑音を同時に除去することができる高い画像復元能力を実現することができるシンプルで実用的な画像復元装置および画像復元方法を提供することである。
【課題を解決するための手段】
【0018】
本発明の画像復元装置は、原画像情報に不必要な情報が混在した劣化画像情報のみから前記原画像情報を推定する画像復元装置であって、時刻nのみの劣化画像情報に対して、時刻nまたは時刻n+1までの情報により前記原画像情報を含む時刻n+1でのシステムの状態量を推定した場合の推定誤差の相関値を算出する相関演算部と、時刻nのみの劣化画像情報に対して、前記相関演算部によって算出された相関値を用いて、時刻n+1までの情報による時刻n+1での前記状態量の最適推定値と、時刻nまでの情報による時刻n+1での前記状態量の最適推定値と、前記劣化画像情報を含む観測量の推定誤差と、の関係を規定するための重み係数を算出する重み係数算出部と、時刻nのみの劣化画像情報に対して、前記重み係数算出部によって算出された重み係数を用いて、時刻nまたは時刻n+1までの情報による時刻n+1での前記状態量の最適推定値を算出する最適推定値算出部と、を有する構成を採る。
【0019】
本発明の画像復元方法は、原画像情報に不必要な情報が混在した劣化画像情報のみから前記原画像情報を推定する画像復元方法であって、時刻nのみの劣化画像情報に対して、時刻nまたは時刻n+1までの情報により前記原画像情報を含む時刻n+1でのシステムの状態量を推定した場合の推定誤差の相関値を算出する相関演算工程と、時刻nのみの劣化画像情報に対して、前記相関演算工程で算出した相関値を用いて、時刻n+1までの情報による時刻n+1での前記状態量の最適推定値と、時刻nまでの情報による時刻n+1での前記状態量の最適推定値と、前記劣化画像情報を含む観測量の推定誤差と、の関係を規定するための重み係数を算出する重み係数算出工程と、時刻nのみの劣化画像情報に対して、前記重み係数算出工程で算出した重み係数を用いて、時刻nまたは時刻n+1までの情報による時刻n+1での前記状態量の最適推定値を算出する最適推定値算出工程と、を有するようにした。
【発明の効果】
【0020】
本発明によれば、ぼけと雑音を同時に除去することができる高い画像復元能力を実現することができるシンプルで実用的な画像復元装置および画像復元方法を得ることができる。
【図面の簡単な説明】
【0021】
【図1】本発明の実施の形態1に係る画像復元装置の構成を示すブロック図
【図2】(A)は、図1の画像復元処理部の構成を示すブロック図、(B)は、図2(A)の第1復元処理部の構成を示すブロック図
【図3】画像の劣化の過程を説明するための図
【図4】発明手法1の状態空間モデルのシステム構成を表すブロック線図
【図5】発明手法1の状態方程式の定式化の具体例を説明するための図であり、特に、(A)は、状態空間モデルの状態過程を示す図、(B)は、処理対象ブロックおよびその時間変化の一例を示す図、(C)は、状態方程式の具体例を示す図
【図6】発明手法1の状態方程式の構成を説明するための図
【図7】発明手法1の観測方程式の定式化の具体例を説明するための図であり、特に、(A)は、状態空間モデルの観測過程を示す図、(B)は、観測方程式の具体例を示す図
【図8】従来一般的な観測方程式の構成(一般例)を説明するための図であり、特に、(A)は、従来の観測方程式の定義を示す図、(B)は、観測量に影響する状態量の範囲を視覚的に示す図
【図9】発明手法1の観測方程式の構成を説明するための図であり、特に、(A)は、発明手法1の観測方程式の定義を示す図、(B)は、観測量に影響する状態量の範囲を視覚的に示す図
【図10】観測遷移行列を構成する係数の割り当て方法を説明するための図
【図11】図9(A)に示す観測遷移行列に基づく、各観測量に影響する状態量の範囲を、割り当てられた係数と共に示す図
【図12】発明手法1のアルゴリズムの一例を示す図
【図13】図12のアルゴリズムを実行する処理手順を示すフローチャート
【図14】発明手法1を視覚的にまとめた説明図
【図15】シミュレーション条件を説明するための図
【図16】原画像「カメラマン」を説明するための図
【図17】従来手法との対比において原画像「カメラマン」に対する発明手法1のシミュレーション結果(視覚評価)を示す図
【図18】図17中の破線の円の周辺領域を拡大した図
【図19】図17と同様のシミュレーション結果(視覚評価)を示す図
【図20】図19中の破線の円の周辺領域を拡大した図
【図21】従来手法との対比において原画像「レナ」に対する発明手法1のシミュレーション結果(視覚評価)を示す図
【図22】図21中の破線の円の周辺領域を拡大した図
【図23】従来手法との対比において原画像に対する発明手法1のシミュレーション結果(客観評価)を示す図
【図24】従来手法との対比において原画像に対する発明手法1のシミュレーション結果(主観評価)を示す図
【図25】図1の画像復元装置の動作の一例を示すフローチャート
【図26】図25のステップS2300の検出方法の一例を説明するための図
【図27】復元モードの実装の一例を説明するための図
【図28】復元モードの実装の他の例を説明するための図
【図29】劣化画像モデルを視覚的に説明するための図
【図30】発明手法1の状態空間モデルを視覚的に説明するための図
【図31】発明手法1の状態過程を視覚的に説明するための図
【図32】発明手法1の観測過程を視覚的に説明するための図
【図33】発明手法1の状態方程式と観測方程式をまとめて示す説明図
【図34】発明手法1のアルゴリズムを説明するための図
【図35】図34のアルゴリズムを実行する処理手順を示すフローチャート
【図36】発明手法1の課題を視覚的に説明するための図
【図37】発明手法2の状態方程式の構成を視覚的に示す説明図であり、特に、(A)は、発明手法1の状態方程式の具体例を示す図、(B)は、時刻の変化に伴う処理対象ブロックの変化を示す図
【図38】発明手法2の観測方程式の構成を視覚的に示す説明図
【図39】発明手法2の観測方程式の構成を視覚的に示す他の説明図
【図40】発明手法2の状態方程式と観測方程式をまとめて示す説明図
【図41】発明手法2のアルゴリズムを説明するための図
【図42】図41のアルゴリズムを実行する処理手順を示すフローチャート
【図43】発明手法1と発明手法2のアルゴリズムを比較した説明図
【図44】発明手法2の演算量評価を説明するための図
【図45】発明手法2の復元処理速度評価を説明するための図
【図46】シミュレーション条件を説明するための図
【図47】比較対象手法との対比において原画像「カメラマン」に対する発明手法2のシミュレーション結果(視覚評価)を示す図
【図48】図47中の破線の円の周辺領域を拡大した図
【図49】比較対象手法との対比において原画像「レナ」に対する発明手法2のシミュレーション結果(視覚評価)を示す図
【図50】図49中の破線の円の周辺領域を拡大した図
【図51】比較対象手法との対比において原画像に対する発明手法2のシミュレーション結果(客観評価)を示す図
【図52】発明手法2の課題を視覚的に説明するための図
【図53】発明手法3の観測方程式の構成を視覚的に示す説明図
【図54】発明手法3の状態空間モデルを視覚的に説明するための図
【図55】発明手法3の状態方程式と観測方程式をまとめて示す説明図
【図56】発明手法3のアルゴリズムを説明するための図
【図57】図56のアルゴリズムを実行する処理手順を示すフローチャート
【図58】発明手法3の演算量評価を説明するための図
【図59】発明手法3の復元処理速度評価を説明するための図
【図60】比較対象手法との対比において原画像「カメラマン」に対する発明手法3のシミュレーション結果(視覚評価)を示す図
【図61】図60中の破線の円の周辺領域を拡大した図
【図62】比較対象手法との対比において原画像「レナ」に対する発明手法3のシミュレーション結果(視覚評価)を示す図
【図63】図62中の破線の円の周辺領域を拡大した図
【図64】比較対象手法との対比において原画像に対する発明手法3のシミュレーション結果(客観評価)を示す図
【図65】1次元の画像復元手法の特徴を説明するための図
【図66】2次元の画像復元手法の特徴を説明するための図
【図67】発明手法4の状態空間モデルを視覚的に説明するための図
【図68】発明手法4の状態過程を視覚的に説明するための図
【図69】発明手法4の状態方程式の構成を視覚的に示す説明図
【図70】発明手法4の観測過程を視覚的に説明するための図
【図71】発明手法4の観測方程式の構成を視覚的に示す説明図
【図72】発明手法4の観測遷移行列の具体例を示す図
【図73】図72の観測遷移行列の構成を視覚的に示す説明図
【図74】図72の観測遷移行列を構成する係数の割り当て方法の具体例を示す図
【図75】図72の観測遷移行列を構成する係数の割り当て方法の具体例を示す他の図
【図76】図72の観測遷移行列を構成する係数の割り当て方法の具体例を示すさらに他の図
【図77】発明手法4の状態方程式と観測方程式をまとめて示す説明図
【図78】発明手法4のアルゴリズムを説明するための図
【図79】図78のアルゴリズムを実行する処理手順を示すフローチャート
【図80】発明手法4の演算量評価を説明するための図
【図81】発明手法4の復元処理速度評価を説明するための図
【図82】比較対象手法との対比において原画像「カメラマン」に対する発明手法4のシミュレーション結果(視覚評価)を示す図
【図83】図82中の破線の円の周辺領域を拡大した図
【図84】比較対象手法との対比において原画像「レナ」に対する発明手法4のシミュレーション結果(視覚評価)を示す図
【図85】図84中の破線の円の周辺領域を拡大した図
【図86】比較対象手法との対比において原画像に対する発明手法4のシミュレーション結果(客観評価)を示す図
【図87】発明手法4の課題を視覚的に説明するための図
【図88】発明手法5の状態方程式の構成を視覚的に示す説明図であり、特に、(A)は、発明手法4の状態方程式の具体例を示す図、(B)は、時刻の変化に伴う処理対象ブロックの変化を示す図
【図89】発明手法5の観測方程式の構成を視覚的に示す説明図
【図90】発明手法5の観測方程式の構成を視覚的に示す他の説明図
【図91】発明手法5の状態方程式と観測方程式をまとめて示す説明図
【図92】発明手法5のアルゴリズムを説明するための図
【図93】図92のアルゴリズムを実行する処理手順を示すフローチャート
【図94】発明手法4と発明手法5のアルゴリズムを比較した説明図
【図95】発明手法5の演算量評価を説明するための図
【図96】発明手法5の復元処理速度評価を説明するための図
【図97】比較対象手法との対比において原画像「カメラマン」に対する発明手法5のシミュレーション結果(視覚評価)を示す図
【図98】図97中の破線の円の周辺領域を拡大した図
【図99】比較対象手法との対比において原画像「レナ」に対する発明手法5のシミュレーション結果(視覚評価)を示す図
【図100】図99中の破線の円の周辺領域を拡大した図
【図101】比較対象手法との対比において原画像に対する発明手法5のシミュレーション結果(客観評価)を示す図
【図102】発明手法5の課題を視覚的に説明するための図
【図103】発明手法6の観測方程式の構成を視覚的に示す説明図
【図104】発明手法6の状態空間モデルを視覚的に説明するための図
【図105】発明手法6の状態方程式と観測方程式をまとめて示す説明図
【図106】発明手法6のアルゴリズムを説明するための図
【図107】図106のアルゴリズムを実行する処理手順を示すフローチャート
【図108】発明手法6の演算量評価を説明するための図
【図109】発明手法6の復元処理速度評価を説明するための図
【図110】比較対象手法との対比において原画像「カメラマン」に対する発明手法6のシミュレーション結果(視覚評価)を示す図
【図111】図110中の破線の円の周辺領域を拡大した図
【図112】比較対象手法との対比において原画像「レナ」に対する発明手法6のシミュレーション結果(視覚評価)を示す図
【図113】図112中の破線の円の周辺領域を拡大した図
【図114】比較対象手法との対比において原画像に対する発明手法6のシミュレーション結果(客観評価)を示す図
【図115】本発明の実施の形態2に係る画像復元装置の構成を示すブロック図
【図116】図115の復元モード指定部によって指定可能な復元モードの具体例を示す図
【図117】図115の画像復元装置の動作の一例を示すフローチャート
【発明を実施するための形態】
【0022】
以下、本発明の実施の形態について、図面を参照して詳細に説明する。
【0023】
(実施の形態1)
図1は、本発明の実施の形態1に係る画像復元装置の構成を示すブロック図である。ここでは、本発明の画像復元方法を適用した好適な画像復元装置について、特に各種の用途に対応可能な汎用的な画像復元装置を例にとって説明する。
【0024】
図1に示す画像復元装置100は、コンピュータで構成されており、大別して、画像入力装置110、入力インタフェース部120、操作部130、内部インタフェース部140、記憶部150、画像復元処理部160、出力インタフェース部170、および画像出力装置180を有する。
【0025】
画像入力装置110は、復元処理の対象となる画像データ(劣化画像)をデジタルデータとしてコンピュータに入力するための入力装置である。入力する画像は、静止画像でも動画像でもよい。画像入力装置110としては、例えば、カメラ112やスキャナ114、記録メディア116、モデム118などを用いることができる。カメラ112は、撮像機能を有するすべての装置を意味し、例えば、デジタルカメラ(デジタルスチルカメラとデジタルビデオカメラ)のほかに、カメラ機能を搭載した携帯電話や、防犯カメラ(監視カメラ)、画像診断を行うための医療機器(内視鏡、レントゲン、エコー、CT、MRIなど)などを含むことができる。スキャナ114は、代表的な画像入力装置の1つであり、ネガやポジといったファイルからの読み込みを専用に行うフィルムスキャナなども含んでいる。記録メディア116は、画像データを記録可能な記録媒体を広く意味し、例えば、磁気ディスク(HDDやFDなど)や、光ディスク(CDやDVD、BDなど)、光磁気ディスク(MO)、フラッシュメモリ(メモリカードやUSBメモリなど)などである。モデム118は、外部の通信ネットワーク(例えば、電話回線やLAN、インターネットなど)と接続するための装置である。なお、画像入力装置110の種類は、画像復元装置100の用途に応じて適宜選択すればよい。
【0026】
入力インタフェース部120は、画像入力装置110から提供される画像データをコンピュータで処理可能なデータ形式に変換するなどの入力処理を行う。図示しないが、入力インタフェース部120は、画像入力装置110の種類に応じて別個独立に設けられる。例えば、記録メディア116の入力インタフェース部120は、ドライブとも呼ばれ、記録メディア116の種類に応じて様々な種類のドライブが使用可能である。なお、ドライブは、記録メディアを読み書きする装置であり、記録メディアに関する限り、通常は、入力インタフェース部120と出力インタフェース部170は一体化されている。また、通常、モデム118は、画像入力装置110としても画像出力装置180としても機能しうるため、モデム118に関しても、通常、入力インタフェース部120と出力インタフェース部170は一体化されている。入力インタフェース部120は、コンピュータ本体の内部に格納された内蔵カード(ボード)であってもよいし、内部インタフェース部140を介して接続された外部設置型機器であってもよい。
【0027】
なお、画像入力装置110が画像情報をアナログデータとして出力する場合、対応する入力インタフェース部120は、サンプリング部およびA/D変換部(共に図示せず)を有する。サンプリング部は、所定のサンプリング周波数で、入力されたアナログ信号をサンプリング処理し、A/D変換部に出力する。サンプリング周波数は、復元処理対象(情報源)の種類に応じて変更可能である。A/D変換部は、サンプリングされた信号の振幅値を所定の分解能でA/D変換処理する。
【0028】
操作部130は、典型的にはキーボートやマウス、タッチパネルなどであるが、音声認識装置などを用いてもよい。使用者は、操作部130を用い、例えば、後述するディスプレイ182で確認をしながらコンピュータを操作することができる。また、操作部130は、パラメータ設定部132、領域指定部134、および復元モード指定部136を有する。パラメータ設定部132は、使用者の入力操作により、後で詳述する本実施の形態における画像復元処理に必要な各種パラメータの値を設定する。領域指定部134は、使用者の入力操作により、入力された画像に対して画像復元処理の対象となる領域(画像の特定の範囲)を指定する。復元モード指定部136は、使用者の入力操作により、後述する復元モードを指定する。
【0029】
内部インタフェース部140は、コンピュータ本体の内部にあって、入力インタフェース部120、操作部130、記憶部150、画像復元処理部160、および出力インタフェース部170を相互に接続する機能を有する。コンピュータの内部における各種信号のやり取りは、内部インタフェース部140を介して行われる。
【0030】
記憶部150は、主記憶装置152および補助記憶装置154を有する。主記憶装置152は、コンピュータの本体を構成する一要素であって、主にプログラムとデータを格納する。補助記憶装置154は、主記憶装置152の容量の不足を補う記憶装置である。補助記憶装置154は、典型的にはハードディスク(HD)であるが、CD−ROMやDVD、SSD(Solid State Drive)、フラッシュメモリなどの可搬性のあるものであってもよいし、また、それらの組み合わせであってもよい。本実施の形態における画像復元処理を実行させるプログラム(画像復元アルゴリズム)は、予め記憶部150(主記憶装置152または補助記憶装置154)に格納されていてもよいし、記録メディア116から両インタフェース部120、140を介して記憶部150にインストールしたり、モデム118および両インタフェース部120、140を介して外部から記憶部150にダウンロードしたりしてもよい。
【0031】
なお、画像入力装置110から画像データを取り込み、取り込んだ画像データに対して画像復元処理を行い、復元処理後の画像データを画像出力装置180から取り出すという一連の処理を実行するためには、取り込んだ画像データを格納しておく記憶領域と、データ処理の上で一時的に必要な記憶領域(作業領域または作業用メモリと呼ばれる)と、出力すべき画像データを格納しておく記憶領域とが必要である。これらの記憶領域は、主記憶装置152または補助記憶装置154の上に取られるが、ここでは、説明の便宜上、復元処理後の画像データを後述するディスプレイ182に出力する場合を想定して、特に表示メモリ156を図面上に別途表記している。
【0032】
画像復元処理部160は、本発明の特徴的な構成要素であり、後述する内蔵された画像復元アルゴリズムを実行する。例えば、カルマンフィルタを用いた従来の画像復元手法では、2段階処理、つまり、AR次数を決定し、AR係数を推定した後、この推定したAR係数を用いて状態空間モデル(状態方程式と観測方程式)を構成してカルマンフィルタを実行することによって、画像復元を実現しているのに対し、本発明の画像復元方法(以下「発明手法」という)では、状態方程式と観測方程式から構成される新しい予測法によって画像復元を実現している。具体的には、発明手法では、ARシステムのコンセプトを必要としない新しい状態空間モデル(状態方程式と観測方程式)を構成することによって、より具体的には、状態方程式をクリアな画像情報(原画像情報)のみで構成し、かつ、観測方程式を劣化画像情報、クリアな画像情報(原画像情報)、ぼけ情報、および雑音で構成した新たな状態空間モデルを用いることによって、画像復元を実現している。画像復元の対象となる画像は、静止画像でも動画像でもよい。なお、本発明の画像復元方法については、後で詳細に説明する。
【0033】
画像出力装置180は、コンピュータ(画像復元処理部160)で復元処理された画像データ(復元画像)を所定の出力形態で外部に取り出すための出力装置である。画像出力装置180としては、例えば、ディスプレイ182やプリンタ184、記録メディア186、モデム188などを用いることができる。記録メディア186およびモデム188は、画像入力装置110としての記録メディア116およびモデム118とそれぞれ共用してもよい。なお、画像出力装置180の種類は、画像復元装置100の用途に応じて適宜選択すればよい。
【0034】
出力インタフェース部170は、コンピュータ(画像復元処理部160)で復元処理された画像データ(復元画像)を画像出力装置180に出力可能なデータ形式に変換するなどの出力処理を行う。図示しないが、出力インタフェース部170は、画像出力装置180の種類に応じて別個独立に設けられる。上記のように、記録メディアおよびモデムに関しては、通常、入力インタフェース部120と出力インタフェース部170は一体化されている。出力インタフェース部170も、入力インタフェース部120と同様に、コンピュータ本体の内部に格納された内蔵カード(ボード)であってもよいし、内部インタフェース部140を介して接続された外部設置型機器であってもよい。
【0035】
図2(A)は、図1の画像復元処理部160の構成を示すブロック図であり、図2(B)は、図2(A)の第1復元処理部160aの構成を示すブロック図である。
【0036】
本実施の形態では、図2(A)に示すように、画像復元処理部160は、第1復元処理部160aと第2復元処理部160bを有する。第1復元処理部160aは、本発明の画像復元方法を実施する。第2復元処理部160bは、本発明の画像復元方法以外の画像復元方法、例えば、ウィナーフィルタを用いた画像復元手法やカルマンフィルタを用いた画像復元手法などの従来の画像復元方法を実施する。これによって、第1復元処理部160aと第2復元処理部160bとを用いた画像復元処理が可能となる。なお、以下では、第2復元処理部160bで実施する従来の画像復元方法として、例えば、ウィナーフィルタを用いた画像復元手法(以下「従来手法1」という)を例にとって説明する。
【0037】
第1復元処理部160aは、図2(B)に示すように、初期設定部162、相関演算部164、重み係数算出部166、および最適推定値算出部168を有する。第1復元処理部160aは、各部162〜168の協働の下、内蔵された画像復元アルゴリズム(発明手法)を実行することによって、取り込んだ画像データ(劣化画像)のみからクリアな画像(原画像)を復元する。このとき、初期設定部162は、発明手法の画像復元アルゴリズムの初期設定を行い、相関演算部164は、原画像(所望情報、クリアな画像)の推定誤差の相関演算を行い、重み係数算出部166は、原画像(所望情報)の最適推定値の算出に必要な重み係数の算出を行い、最適推定値算出部168は、原画像(所望情報)の最適推定値の算出を行う。なお、各部162〜168の具体的な処理内容は、後で詳細に説明する。
【0038】
以下、第1復元処理部160aで行われる画像復元処理動作について詳細に説明するが、ここでは、まず、画像復元の理論となる劣化画像モデルについて説明しておく。
【0039】
劣化した画像は、一般的に、原画像とぼけの点拡がり関数(PSF:Point Spread Function、以下単に「ぼけ関数」ともいう)の畳み込み(コンボリューション)に雑音を加えたモデルによって得られる。すなわち、イメージングシステム(対象物を撮像する装置や画像を生成するシステムなど)によって得られる劣化画像g(x,y)は、対象物を表す原画像をf(n,m)、ぼけの点拡がり関数をh(x,y)、雑音をn(x,y)でそれぞれ表すと、次の式(1)で表される。なお、ぼけの点拡がり関数h(x,y)は、撮影条件なども含めたイメージングシステムの特性を表している。
【0040】
【数1】

【0041】
図3は、画像の劣化の過程を説明するための図である。
【0042】
例えば、図3の左端に示す原画像にぼけが生じると、図3の中央に示す画像となり、この画像にさらに雑音が加わると、図3の右端に示す劣化画像となる。ぼけは、ある画素がその周囲の画素に影響することによって生じるが、雑音は、画素に関係なく発生する。上記のように、カメラなどで撮影された画像の場合、ぼけは、手ぶれや焦点ずれなどに起因して生じ、雑音は、暗電流や熱雑音などに起因して不可避的に発生する。図3に示すように、劣化画像の復元は、劣化画像からぼけと雑音を取り除くことにほかならない。
【0043】
発明手法では、上記のように、ARシステムのコンセプトを用いないようにするために、新しい状態空間モデル(状態方程式と観測方程式)を構成する。すなわち、クリアな画像情報(原画像情報)のみを用いて状態方程式を構成し、かつ、劣化画像情報、クリアな画像情報(原画像情報)、ぼけ情報、および雑音を用いて観測方程式を構成する。なお、本実施の形態では、発明手法として、6つの具体的な手法を提示する。ここでは、便宜上、第1の発明手法から第6の発明手法を、発明手法1、発明手法2、発明手法3、発明手法4、発明手法5、発明手法6とそれぞれ称し、順次説明する。
【0044】
発明手法1〜6は、いずれも、AR係数を用いない画像復元手法であって、ぼけと雑音を同時に除去することができる、シンプルで実用的な高性能の画像復元手法である。ただし、発明手法1〜6は、互いに、アルゴリズムの違いに応じて、演算量や復元処理速度、画像復元能力などが異なる。ここでは、各発明手法1〜6を詳細に説明する前に、各発明手法1〜6の特徴、長所、短所をまとめておく。
【0045】
発明手法1〜3は、1次元の画像復元手法であり、発明手法4〜6は、発明手法1〜3をそれぞれ2次元に拡張した2次元の画像復元手法である。発明手法1は、本発明に係る画像復元方法の基本型である。発明手法1は、AR係数を用いない画像復元手法であって、有色性駆動源を持つ状態空間モデルを構成している。発明手法1は、従来技術に対して、高性能な劣化画像復元が可能であるという長所を有する反面、他の発明手法2、3に比べて、演算量が多いという短所を有する。発明手法2は、発明手法1に対して、注目画素領域のみ考慮した状態空間モデルを再構成した手法であって、注目画素領域のみ用いることで、演算量を大幅に軽減している。発明手法2は、高速(低演算量)な劣化画像復元が可能であるという長所を有する(画像復元能力は発明手法1とほぼ同等である)。発明手法3は、発明手法1、2に対して、有色駆動源の影響も考慮した状態空間モデルを再構成した手法である。発明手法3は、演算量の軽減、復元処理速度の向上、有色性駆動源にも適用可能、および画像復元能力の向上といった長所を有し、画像復元能力を向上しながら高速(低演算量)な劣化画像復元が可能な手法である。
【0046】
上記のように、発明手法1は、AR係数を用いない1次元の画像復元手法であって、従来技術に対して、高性能な劣化画像復元が可能であるという特徴を有する。発明手法4〜6は、2次元の画像復元手法であって、発明手法1に対して、柔軟な処理効率の向上が図られるという特徴を有する。発明手法4は、発明手法1を2次元に拡張した画像復元手法であって、注目画素領域と周辺画素領域をそれぞれ任意にとることにより処理効率の向上が図られるという長所を有する反面、他の発明手法1〜3、5、6に比べて、演算量が多いという短所を有する。発明手法5は、発明手法4に対して、注目画素領域のみ考慮した状態空間モデルを再構成した手法であって、注目画素領域のみ用いることで、発明手法4に対して、演算量を大幅に軽減している。発明手法5は、高速(低演算量)な劣化画像復元が可能であるという長所を有する反面、発明手法4よりも若干復元精度が低下するという短所を有する。発明手法6は、発明手法4、5に対して、有色駆動源の影響も考慮した状態空間モデルを再構成した手法であって、発明手法5の短所を改善したものである。発明手法6は、良好な画像復元能力、有色性駆動源にも適用可能、演算量の軽減、復元する画素に適した復元処理が可能といった長所を有し、処理効率性を考慮した実用的な劣化画像復元手法である。特に、発明手法1〜3を2次元に拡張した発明手法4〜6は、注目画素領域と周辺画素領域をそれぞれ任意に選択することができるため、劣化画像の状態と使用者の要望(処理速度または処理品質)に応じて、1次元の発明手法1〜3に適応することも可能であるし、あるいは、若干の変更により2次元の発明手法4〜6に適応することも可能である、という特徴を有する。
【0047】
具体的には、発明手法では、新しい状態空間モデル(状態方程式と観測方程式)を構成し、この新しい状態空間モデルは、次の式(2)で表される。ただし、状態方程式において、ベクトルxp1(n)は状態ベクトル(原画像情報)、行列Φp1は状態遷移行列、ベクトルδp1(n)は駆動源ベクトルである。また、観測方程式において、ベクトルyp1(n)は観測ベクトル(観測画像情報、劣化画像情報)、行列Mp1は観測遷移行列、ベクトルεp1(n)は観測雑音ベクトルである。また、「n」とは、装置の時刻nである。時刻nは、注目画素(または注目画素領域)と複数の周辺画素(または周辺画素領域)とからなる処理対象ブロックの処理の順番(n番目)を表している。なお、ここでは、説明の便宜上、発明手法として、例えば、基本型である発明手法1を提示する。以下の説明において、添え字「p1」は、発明手法1に係るものであることを示している。また、ここでは、注目画素領域および周辺画素領域における「画素領域」は、1つの画素のみを含む場合と、複数の画素を含む場合の両方を意味する。すなわち、ぼけと雑音が含まれる劣化画像を復元する際の復元処理対象として、1つの画素を「注目画素」と称し、少なくとも1つ以上の画素を含む領域を「注目画素領域」と称する。また、注目画素に隣接する画素を「周辺画素」と称し、注目画素に隣接する画素を含む領域、または、注目画素に隣接する画素および注目画素に隣接しないがその周辺の画素を含む領域を、「周辺画素領域」と称する。
【0048】
【数2】

【0049】
図4は、この状態空間モデルにおけるシステム構成をブロック線図で表したものである。
【0050】
図4に示すように、この状態空間モデルは、状態過程と観測過程とで構成される。状態過程は、状態方程式によって記述され、観測過程は、観測方程式によって記述される。図4において、「201」は時刻nにおける状態ベクトルxp1(n)、「202」は時刻n+1における状態ベクトルxp1(n+1)、「203」は時刻nにおける観測ベクトルyp1(n)、「204」は時刻nにおける観測雑音ベクトルεp1(n)、「205」は時刻n+1における駆動源ベクトルδp1(n+1)、「206」は状態遷移行列Φp1、「207」は観測遷移行列Mp1である。式(2)における状態方程式は、観測対象のシステムを状態空間モデルで記述したものであり、内部状態つまり状態変数(ここでは、状態ベクトルxp1(n)、xp1(n+1))の時間に対する生成過程を表している。また、式(2)における観測方程式は、何らかの観測装置を通じて観測する過程を記述したものであり、観測結果(ここでは、観測ベクトルyp1(n))が、被観測量つまり入力(ここでは、状態ベクトルxp1(n))に依存して生成される様子を示している。
【0051】
図5は、発明手法1の状態方程式の定式化の具体例を説明するための図であり、特に、図5(A)は、状態空間モデルの状態過程を示す図、図5(B)は、処理対象ブロックおよびその時間変化の一例を示す図、図5(C)は、状態方程式の具体例を示す図である。
【0052】
発明手法1では、劣化画像の復元処理は、復元対象である注目画素のみを用いるのではなく、その周辺の画素も含めて処理を行う。すなわち、J×Jの大きさの画像のある注目画素を中心にして、その周辺K×K(ただし、J>K)の領域を処理対象ブロックとし、この処理対象ブロックの中の画素情報を用いて、画像の復元処理を行う。このように、処理対象ブロックとは、画像復元においてK×Kの画素を用いて処理を行う範囲を意味する。
【0053】
例えば、J=256、K=3とした場合、図5に示すように、256×256の大きさの画像のある注目画素を中心にして、その周辺3×3の画素領域が処理対象ブロックとなる。図面上、処理対象ブロックは、網かけされている。画像復元処理は、この3×3の処理対象ブロックの中の画素情報を用いて行われる。単に説明の便宜上、図5(B)に示すように、256×256の大きさの画像の一部の画素に1〜36の番号を付したとすると、時刻nにおける状態ベクトルxp1(n)に対応する処理対象ブロックは、番号が1、2、3、7、8、9、13、14、15の3×3の9個の画素を構成要素とし、次の時刻n+1における状態ベクトルxp1(n+1)に対応する処理対象ブロックは、番号が7、8、9、13、14、15、19、20、21の3×3の9個の画素を構成要素としている。
【0054】
このとき、式(2)における状態方程式は、図5(C)にも示す次の式(3)で定義される。ここで、状態ベクトルxp1(n)は、原画像情報を含む状態量として3×3の処理対象ブロックに含まれる9個の画素情報x(n)、x(n)、x(n)、x(n)、x(n)、x(n)、x(n)、x(n)、x(n)を要素とする9×1次のベクトルである。また、状態遷移行列Φp1は、式(3)で定義された9×9次の行列であり、駆動源ベクトルδp1(n)は、同じく式(3)で定義された9×1次のベクトルである。
【0055】
【数3】

【0056】
図6は、式(3)で表される状態方程式の構成を説明するための図である。
【0057】
式(3)で表される状態方程式の特徴は、状態遷移行列Φp1の一部の要素を「1」に設定し、残りの要素をすべて「0」に設定した点、および、駆動源ベクトルδp1(n)の一部の要素を有色信号である状態量x(n)(i=7,8,9)で表した点にある。これは、図6に示すように、現在の時刻n+1における処理対象ブロックと過去の時刻nにおける処理対象ブロックとが重なる部分において、現在の状態量(原画像の画素情報)x(n+1)、x(n+1)、x(n+1)、x(n+1)、x(n+1)、x(n+1)と過去の状態量(原画像の画素情報)x(n)、x(n)、x(n)、x(n)、x(n)、x(n)とを対応させるためである。この結果、xp1(n)とxp1(n+1)の関係を表す状態方程式は、クリアな画像からなる状態ベクトルxp1(n)と、0と1からなる状態遷移行列Φp1と、有色信号であるクリアな画像からなる駆動源ベクトルδp1(n)とから構成されるため、所望の状態量(原画像の画素情報)、つまり、クリアな画像情報(原画像情報)のみで構成されることになる。
【0058】
図7は、発明手法1の観測方程式の定式化の具体例を説明するための図であり、特に、図7(A)は、状態空間モデルの観測過程を示す図、図7(B)は、観測方程式の具体例を示す図である。
【0059】
図5および図6の例に対応して、式(2)における観測方程式は、図7(B)にも示す次の式(4)で定義される。ここで、観測ベクトルyp1(n)は、劣化画像情報を含む観測量として3×3の処理対象ブロックに含まれる9個の画素情報y(n)、y(n)、y(n)、y(n)、y(n)、y(n)、y(n)、y(n)、y(n)を要素とする9×1次のベクトルである。また、観測遷移行列Mp1は、式(4)で定義された9×9次の行列であって、画像の劣化モデルにおけるぼけの点拡がり関数(PSF)に対応している。観測遷移行列Mp1を構成する各要素hr,s(r,sはhの座標であり、r,s=−1,0,1)は、既知であって、予めデータ化して適切に設定されている。また、観測雑音ベクトルεp1(n)は、3×3の処理対象ブロックに含まれる9個の画素に対応する観測雑音ε(n)、ε(n)、ε(n)、ε(n)、ε(n)、ε(n)、ε(n)、ε(n)、ε(n)を要素とする9×1次のベクトルである。
【0060】
【数4】

【0061】
図8は、従来一般的な観測方程式の構成(一般例)を説明するための図であり、図9は、式(4)で表される発明手法1の観測方程式の構成を説明するための図である。特に、図8(A)は、従来一般的な観測方程式の定義を示す図、図8(B)は、観測量に影響する状態量の範囲を視覚的に示す図である。また、図9(A)は、発明手法1の観測方程式の定義を示す図、図9(B)は、観測量に影響する状態量の範囲を視覚的に示す図である。なお、図8および図9の説明においては、便宜上、観測雑音ベクトルεp1(n)は省略する。
【0062】
画像においてぼけに基づく劣化とは、上記のように、ある画素がその周囲の画素の影響を受けて劣化することである。そこで、従来は、一般的に、図8(A)に示すように観測方程式を定義していた。この場合、ぼけは、周囲のすべての画素の影響を受けることに起因して生じる。すなわち、処理対象ブロックの中の各画素は、処理対象ブロックの中のすべての画素の影響を受けることになる。言い換えれば、各観測量y(i=1,2,…,9)は、すべての状態量x(i=1,2,…,9)の影響を受けて成り立つ。これを視覚的に示すと、例えば、図8(B)に示すようになる。すなわち、図8(B)に示すように、例えば、観測量がyの場合もyの場合も、3×3(=9)個のすべての状態量x〜xの影響を受けることになる。
【0063】
これに対し、式(4)で表される発明手法1の観測方程式の特徴は、観測遷移行列Mp1の一部の要素を、3×3(=9)個の係数hr,s(座標r,s=−1,0,1)を用いて規則的に設定し、残りの要素をすべて「0」に設定した点にある。例えば、図9(A)は、観測方程式を展開して得られる一部の式(y(n)、y(n)に関する式)を示している。観測遷移行列Mp1に対する係数hr,sの割り当て方法は、次の通りである。すなわち、9個の係数hr,s(r,s=−1,0,1)を図10に示すようにマトリックス状に配列し(以下「係数マトリックス」という)、係数マトリックスの中央(つまり、係数h0,0)の位置を観測量y(i=1,2,…,9)の注目画素iの位置に合わせたときの、その係数マトリックスに従って、観測遷移行列Mp1の一部の要素に係数hr,s(r,s=−1,0,1)を割り当てる。これを視覚的に示すと、例えば、図9(B)に示すようになる。すなわち、図9(B)に示すように、例えば、観測量がyの場合は、3×3(=9)個の状態量x〜xのうち、5個のx〜xの影響を受け、観測量がyの場合は、6個のx〜xの影響を受けることになる。
【0064】
図11は、図9(A)に示す観測遷移行列Mp1に基づく、各観測量yに影響する状態量xの範囲を、割り当てられた係数hr,s(r,s=−1,0,1)と共に示す図である。なお、ここでは、説明の便宜上、図11(A)〜図11(I)の各図において、左側の処理対象ブロックには、注目する観測量yのみを示し、右側の処理対象ブロックには、影響を与える各状態量xの場所に、この状態量xに乗算されるhr,sのみを示している。
【0065】
式(4)で表される観測方程式では、観測量yは、図11(A)に示すように、5個の状態量x〜xの影響を受け、観測量yは、図11(B)に示すように、6個の状態量x〜xの影響を受け、観測量yは、図11(C)に示すように、7個の状態量x〜xの影響を受ける。また、観測量yは、図11(D)に示すように、8個の状態量x〜xの影響を受け、観測量yは、図11(E)に示すように、9個すべての状態量x〜xの影響を受け、観測量yは、図11(F)に示すように、8個の状態量x〜xの影響を受ける。また、観測量yは、図11(G)に示すように、7個の状態量x〜xの影響を受け、観測量yは、図11(H)に示すように、6個の状態量x〜xの影響を受け、観測量yは、図11(I)に示すように、5個の状態量x〜xの影響を受ける。
【0066】
図12は、発明手法1のアルゴリズムの一例を示す図である。なお、このアルゴリズムは、画像の種類を問わず、静止画像のみならず動画像にも適用可能である。
【0067】
発明手法1のアルゴリズムは、図12に示すように、初期設定(Initialization)の過程と反復(Iteration)の過程とに大別され、反復過程は、新しい状態空間モデル(状態方程式と観測方程式)を基に構成されている。反復過程では、1〜6の手順を逐次繰り返す。
【0068】
図13は、図12のアルゴリズムを実行する処理手順を示すフローチャートである。
【0069】
まず、初期設定部162で、初期設定を行う(S1000)。具体的には、初期設定部162において、状態ベクトル、つまり、状態量である所望信号(原画像信号)ベクトルの最適推定値(以下「所望信号の最適推定値ベクトル」という)の初期値xp1(0|0)、所望信号ベクトルの推定誤差(以下「所望信号の推定誤差ベクトル」という)の相関行列の初期値Pp1(0|0)、駆動源ベクトルの共分散の値Rδp1(n)[i,j]、および、観測雑音ベクトルの共分散の値Rεp1(n)[i,j]を、次の式(5)に示すようにそれぞれ設定する。また、図12では図示しないが、状態遷移行列Φp1および観測遷移行列Mp1を、式(3)および式(4)に示すようにそれぞれ設定する。なお、ここでは、時刻nのカウンタの初期値を「0」に設定している。また、ベクトルと行列の要素を示す場合、ベクトルa(n)のi番目の要素をa(n)[i]と表記し、行列A(n)のi行j列の要素をA(n)[i,j]と表記することとする。また、図12において、所望信号の最適推定値ベクトルxp1は、ハットxp1と表記する。
【0070】
【数5】

【0071】
ただし、ベクトル0は、K次元の零ベクトルであり、行列Iは、K次の単位行列である。また、Kは、K×Kの大きさの処理対象ブロック(注目画素領域と周辺画素領域)のその次元である。K=3の場合、K−K=6となり、「i,j>K−K」は、「i,j=7,8,9」に対応している。また、E[δp1(n),δp1(n)]は、駆動源ベクトルδp1(n)の自己相関行列である。E[εp1(n)εp1(n)][i,j]は、観測雑音ベクトルεp1の自己相関行列であり、ここではσε[i,j]と等しく、既知と仮定している。ここでいう「既知」とは、別の任意の方法(アルゴリズム)で求められて与えられることを意味する。もし観測雑音ベクトルεp1(n)が白色雑音であり零平均であれば、σεは所定のサンプル数の分散で与えられる。
【0072】
次に、相関演算部164で、n→(n+1)の誤差共分散行列、つまり、n→(n+1)の推定誤差の相関値(行列)(以下「相関行列」ともいう)を計算する(S1100)。具体的には、相関演算部164において、時刻nまでの情報により時刻n+1での所望信号ベクトルを推定した場合の誤差(所望信号の推定誤差ベクトル)の相関行列Pp1(n+1|n)を計算する。この計算は、ステップS1000でそれぞれ設定した状態遷移行列Φp1および駆動源ベクトルの共分散Rδp1(n+1)[i,j]の値、ならびに、ステップS1000で設定した(n=0の場合)または前回のステップS1600で計算した(n≧1の場合)所望信号の推定誤差ベクトルの相関行列Pp1(n|n)を用いて、次の式(6)により行う。なお、このステップS1100は、図12の反復過程の手順1に相当する。
【0073】
【数6】

【0074】
次に、重み係数算出部166で、重み係数(行列)の計算を行う(S1200)。具体的には、重み係数算出部166において、観測量である観測信号ベクトルの推定誤差(以下「観測信号の推定誤差ベクトル」という)に重み係数(行列)をかけて、時刻nまでの情報による時刻n+1での所望信号の最適推定値ベクトルxp1(n+1|n)を加えたものが、時刻n+1までの情報による時刻n+1での所望信号の最適推定値ベクトルxp1(n+1|n+1)になるような、その重み係数行列Kp1(n+1)を計算する。この計算は、ステップS1100で計算した所望信号の推定誤差ベクトルの相関行列Pp1(n+1|n)、ならびに、ステップS1000でそれぞれ設定した観測遷移行列Mp1および観測雑音ベクトルの共分散Rεp1(n+1)[i,j]の値を用いて、次の式(7)により行う。なお、このステップS1200は、図12の反復過程の手順2に相当する。
【0075】
【数7】

【0076】
次に、最適推定値算出部168で、n→(n+1)の状態量(所望信号)の最適推定値(ベクトル)を計算する(S1300)。具体的には、最適推定値算出部168において、時刻nまでの情報による時刻n+1での所望信号の最適推定値ベクトルxp1(n+1|n)を計算する。この計算は、ステップS1000で設定した状態遷移行列Φp1、および、前回のステップS1400で計算した所望信号の最適推定値ベクトルxp1(n|n)を用いて、次の式(8)により行う。なお、このステップS1300は、図12の反復過程の手順3に相当する。
【0077】
【数8】

【0078】
次に、同じく最適推定値算出部168で、(n+1)→(n+1)の状態量(所望信号)の最適推定値(ベクトル)を計算する(S1400)。具体的には、最適推定値算出部168において、時刻n+1までの情報による時刻n+1での所望信号の最適推定値ベクトルxp1(n+1|n+1)を計算する。この計算は、ステップS1300で計算した所望信号の最適推定値ベクトルxp1(n+1|n)、ステップS1200で計算した重み係数行列Kp1(n+1)、ステップS1000で設定した観測遷移行列Mp1、および時刻n+1における観測信号ベクトルyp1(n+1)を用いて、次の式(9)により行う。なお、このステップS1400は、図12の反復過程の手順4に相当する。
【0079】
【数9】

【0080】
次に、処理を終了するか否かを判断する(S1500)。この判断は、例えば、時刻nが所定のサンプル数Nに達したか否かを判定することによって行う。この判断の結果として、時刻nが所定のサンプル数Nに達していない場合は(S1500:NO)、ステップS1600に進み、時刻nが所定のサンプル数Nに達した場合は(S1500:YES)、ステップS1800に進む。なお、判断の基準は、上記の例に限定されない。例えば、リアルタイムで処理を行う場合は、時刻nが所定のサンプル数Nに達していなくても、サンプルがなくなった時点で処理を終了するようにしてもよい。
【0081】
ステップS1600では、相関演算部164で、(n+1)→(n+1)の誤差共分散行列、つまり、(n+1)→(n+1)の推定誤差の相関値(行列)を計算する。具体的には、相関演算部164において、時刻n+1までの情報により時刻n+1での所望信号ベクトルを推定した場合の誤差(所望信号の推定誤差ベクトル)の相関行列Pp1(n+1|n+1)を計算する。この計算は、ステップS1200で計算した重み係数行列Kp1(n+1)、ステップS1000で設定した観測遷移行列Mp1、およびステップS1100で計算した所望信号の推定誤差ベクトルの相関行列Pp1(n+1|n)を用いて、次の式(10)により行う。なお、このステップS1600は、図12の反復過程の手順5に相当する。
【0082】
【数10】

【0083】
次に、ステップS1700では、時刻nのカウンタを1だけインクリメントして(n=n+1)、ステップS1100に戻る。
【0084】
一方、ステップS1800では、本アルゴリズムの計算結果を出力値として一時保存する。具体的には、ステップS1400で計算した所望信号の最適推定値ベクトルxp1(n+1|n+1)を、本アルゴリズムの出力値として画像復元処理部160内または記憶部150内に一時保存する。
【0085】
図14は、発明手法1を視覚的にまとめた説明図である。このように、発明手法1では、新しい状態空間モデル(状態方程式と観測方程式)を構成しているため、1段階処理で画像の復元処理が可能となる。このことは、他の発明手法2〜6も含めて本発明の大きな特徴の1つである。
【0086】
ここで、カルマンフィルタを用いた従来の画像復元手法との対比において、発明手法1の特徴と効果をまとめて説明しておく。
【0087】
カルマンフィルタを用いた従来の画像復元手法では、上記のように、2段階処理(ステップ1で、AR次数を決定し、AR係数を推定した後、ステップ2で、この推定したAR係数を用いて状態空間モデル(状態方程式と観測方程式)を構成してカルマンフィルタを実行する)によって、画像の復元を実現している。したがって、ステップ1でのAR次数の決定およびAR係数の推定の精度によってステップ2でのカルマンフィルタによる画像復元能力が大きく左右されることが、容易に予想される。これに対して、発明手法1では、ARシステムのコンセプトを必要としない新しい状態空間モデル(状態方程式と観測方程式)を構成し、これを用いた1段階処理の新しい予測法によって、高性能な画像復元を実現している。また、発明手法1では、処理ステップを1段階減らすことができるため、演算量を削減することができ、ひいては、回路規模の縮小やメモリ容量の削減などを図ることができる。
【0088】
また、カルマンフィルタを用いた従来の画像復元手法では、ステップ1でAR係数の推定を行う際に必要となるAR係数の次数の決定が大きな問題点となる。一般にAR係数の次数は状態量に依存しているため、状態量が既知でなければ、理論的にAR係数の次数を正確に決定することは困難である。このことは、状態量が既知でなければならないことを意味しているため、リアルタイム処理は困難となる。また、それにより、正確でないAR係数の次数を用いることになるため、正確なAR係数の推定は困難である。したがって、このことは、カルマンフィルタを用いた従来の画像復元手法の画像復元能力を劣化させる大きな要因となっている。また、何らかの技術でAR次数とAR係数をリアルタイムに推定できるようになったとしても、処理ステップが1段階増えることより演算量の増加は避けられない。これに対して、発明手法1では、ARシステムのコンセプトを必要としないため、このような問題は生じない。
【0089】
また、カルマンフィルタを用いた従来の画像復元手法では、状態量をARシステムを用いて表現することによってモデル化している。このことは、カルマンフィルタを用いた従来の画像復元手法はARシステムでモデル化可能な状態量にしか適用できないことを意味している。すなわち、カルマンフィルタを用いた従来の画像復元手法は、ARシステムでモデル化が困難な状態量には適用することができない。これに対して、発明手法1は、ARシステムのコンセプトを必要としないため、適用対象に関するこのような制約は存在しない。
【0090】
また、カルマンフィルタを用いた従来の画像復元手法では、状態方程式の駆動源が白色信号で、かつ、状態量と観測雑音が無相関であることを仮定した上で、カルマンフィルタ理論を適用している。これに対して、発明手法1では、駆動源が有色信号(クリアな画像)となるものの、状態方程式と観測方程式の特殊な構成によって発明手法1のアルゴリズムを実行することができる。このことは、発明手法1は、カルマンフィルタ理論の一般的な適用条件を考慮しなくても実行可能であることを意味している。すなわち、発明手法1は、カルマンフィルタ理論よりも適用範囲が広いといえる。
【0091】
したがって、ARシステムのコンセプトを必要としない発明手法1は、画像の撮り直しが許されない場面、例えば、医療分野における瞬間的な心臓や肺などの画像復元、汚れなどで劣化した古文書などの復元、および文字・物体認識においても貢献可能な技術である。
【0092】
また、本発明者は、本発明の効果(発明手法1の有効性)を実証するためにシミュレーションを行った。具体的には、発明手法1の画像復元能力を評価するために、(1)視覚評価、(2)客観評価、および(3)主観評価を行った。視覚評価は、原画像と復元画像とを視覚により比較した評価である。客観評価は、数値による評価である。主観評価は、聞き取り調査である。以下、これらを順に説明する。
【0093】
まず、シミュレーションの条件について説明しておく。
【0094】
図15は、シミュレーション条件を説明するための図である。
【0095】
本シミュレーションでは、原画像として、図15に示す2つの画像、つまり、(a)「カメラマン(Cameraman)」と(b)「レナ(Lenna)」を用いた。また、原画像に加えるぼけに対応する点拡がり関数(PSF)のモデルとして、図15に示す2次元ガウス関数を用い、雑音として、加法性白色ガウス雑音を用いた。また、信号対雑音比SNRを20dBとした。また、比較のため、シミュレーションは、同一の上記シミュレーション条件の下で、従来手法1(ウィナーフィルタを用いた画像復元手法)と、従来手法2(射影フィルタを用いた画像復元手法)と、従来手法3(カルマンフィルタを用いた画像復元手法)と、発明手法1に対して行った。例えば、画像「カメラマン」は、図16に示すように、空の部分は画像の分散が変化していないため定常性が強く、人物の頭の部分は画像の分散が変化しているため非定常性が強いといえる。
【0096】
(1)視覚評価
図17は、原画像「カメラマン」に対するシミュレーション結果(視覚評価)を示す図であり、図18は、図17中の破線の円の周辺領域を拡大した図である。
【0097】
図18によく示すように、従来手法1は、ぼけがほとんど取り除かれていない。また、従来手法2は、従来手法1よりもぼけが取り除かれているが、原画像と比べるとまだぼけが残っている。一方、従来手法3は、従来手法1、2に比べて、ぼけは取り除かれてくっきりしているように見えるが、原画像とは異なった復元画像となっている。特に、従来手法3は、原画像と比べて全体的に暗くなっており、空の部分においては、劣化画像よりも劣化してしまっている。
【0098】
これに対して、発明手法1は、図18によく示すように、従来手法1、2、3に比べて原画像を忠実に復元していることが確認できる。すなわち、図18の拡大画像においても、発明手法の有効性を確認することができる。
【0099】
図19は、図17と同様の、原画像「カメラマン」に対するシミュレーション結果(視覚評価)を示す図であり、図20は、図19中の破線の円の周辺領域を拡大した図である。
【0100】
図20によく示すように、カメラの三脚部分に注目すると、従来手法1は、全くぼけが取り除かれていない。また、従来手法2は、従来手法1よりはぼけが取り除かれているが、原画像までは復元されていない。一方、従来手法3は、従来手法1、2に比べるとぼけは取り除かれているように見えるが、画像全体が荒れてしまっており、原画像とかけ離れた復元画像となっていることが確認できる。
【0101】
これに対して、発明手法1は、図20によく示すように、従来手法1、2、3に比べて原画像をより忠実に復元していることが確認できる。すなわち、図20の拡大画像においても、発明手法1の有効性を確認することができる。
【0102】
図21は、原画像「レナ」に対するシミュレーション結果(視覚評価)を示す図であり、図22は、図21中の破線の円の周辺領域を拡大した図である。
【0103】
図22によく示すように、従来手法1は、ほとんどぼけが取り除かれておらず、復元画像が原画像と比べて全体的に白くなっている(輝度が高くなっている)。また、従来手法2は、従来手法1に比べるとぼけは取り除かれているが、原画像まではぼけが取り除かれていない。一方、従来手法3は、従来手法1、2に比べるとぼけが取り除かれておりはっきりとしているように見えるが、劣化画像よりも劣化してしまっており、原画像とはかけ離れた画像になっている。特に肌の部分に注目すると、この点が顕著に表れている。
【0104】
これに対して、発明手法1は、図22によく示すように、従来手法1、2、3に比べて原画像を忠実に復元していることが確認できる。すなわち、図22の拡大画像においても、発明手法1の有効性を確認することができる。
【0105】
(2)客観評価(数値による評価)
図23は、原画像に対するシミュレーション結果(客観評価)を示す図である。
【0106】
ここでは、従来手法1〜3と発明手法1の画像復元能力を数値により評価するため、図23にも示す次の式(11)で表されるPSNR[dB]を用いて画像復元能力を評価した。なお、PSNR(Peak Signal-to-Noise Ratio)は、画像の信号に対する雑音の比であり、数値が大きいほど劣化(ぼけや雑音など)が少なく画像として良好であるといえる。
【0107】
【数11】

【0108】
この場合、図23に示すように、「カメラマン」、「レナ」どちらの画像についても、発明手法1は、従来手法1、2、3よりもSNROUTの数値が大きいことが確認できる。これにより、発明手法1は、客観評価の点でも、従来手法1、2、3に比べて画像復元能力が高いことがわかる。
【0109】
(3)主観評価(聞き取り調査)
図24は、原画像に対するシミュレーション結果(主観評価)を示す図である。
【0110】
ここでは、従来手法1〜3と発明手法1の画像復元能力を評価するために、聞き取り調査による主観評価を行った。画像復元の性能評価は、ACR(絶対範疇評価)に基づく5段階MOS(平均オピニオン値)を用いた聞き取り調査により行った。MOS評価基準は、図24に示す通りである。50人の聴取者が画像復元により得られた画像(図17〜図22参照)を評価した。各々の聴取者は、評価値「1」から「5」を決定する。評価値「5」が最良である。
【0111】
図24に示すように、MOS評価は、「カメラマン」、「レナ」どちらの画像についても、発明手法1は、従来手法1、2、3よりも高い評価を得ていることが確認できる。これにより、発明手法1は、主観評価の点でも、従来手法1、2、3に比べて画像復元能力が高いことがわかる。
【0112】
以上の実験結果により、本発明の画像復元方法(発明手法1)は、従来手法1〜3に比べて高い画像復元能力を発揮していることがわかる。特に、図17〜図22の視覚評価から明らかなように、本発明の画像復元方法は、非定常性が強いエッジ部分などにおいて、従来手法1〜3に比べて著しく復元精度が高いことがわかる。
【0113】
次いで、上記構成を有する画像復元装置100の動作について、図25に示すフローチャートを用いて説明する。なお、図25に示すフローチャートは、記憶部150(主記憶装置152または補助記憶装置154)に制御プログラムとして格納されており、図示しないCPUによって実行される。
【0114】
まず、ステップS2000では、画像入力装置110から復元処理対象の画像データ(劣化画像)を取り込み、記憶部150(主記憶装置152または補助記憶装置154)の所定の記憶領域に格納する。
【0115】
そして、ステップS2050では、ステップS2000で取り込んだ画像データを表示メモリ156に書き込んで、ディスプレイ182に表示する。
【0116】
そして、ステップS2100では、操作部130の領域指定部134で、領域指定処理を行う。具体的には、使用者の入力操作により、ステップS2050でディスプレイ182に表示した画像のうち画像復元処理の対象となる領域(画像の特定の範囲)が指定されると、その指定された領域のデータを生成する。領域の指定は、例えば、画面上、ポインタなどで行われる。なお、使用者により領域の指定が行われなかった場合は、表示された画像全体を指定したものとして取り扱われる。
【0117】
そして、ステップS2150では、拡大処理を行うか否かを判断する。この判断は、例えば、ステップS2100で生成した指定領域データに基づいて行われる。具体的には、指定された領域が表示された画像全体よりも小さい場合には、拡大処理を行うと判断する。この判断の結果として、拡大処理を行う場合は(S2150:YES)、ステップS2200に進み、拡大処理を行わない場合は(S2150:NO)、直ちにステップS2250に進む。
【0118】
ステップS2200では、ステップS2100で指定した領域の拡大処理を行う。具体的には、例えば、指定された領域が表示された画像全体に対応する大きさとなるように拡大処理を行う。拡大処理の結果は、記憶部150(主記憶装置152または補助記憶装置154)の作業用メモリに書き込まれる。なお、拡大処理が行われなかった場合も、その画像データは作業用メモリに書き込まれる。
【0119】
そして、ステップS2250では、時刻nの処理対象ブロック(例えば、3×3の大きさ)を選択する。
【0120】
そして、ステップS2300では、ステップS2250で選択した処理対象ブロックに対して、画像の急激な変化領域が検出されたか否かを判断する。画像の急激な変化領域は、例えば、画像のエッジ部分などに相当する。具体的には、例えば、作業用メモリに書き込まれた画像データ(拡大処理の有無を問わず)に対して、ステップS2250で選択した処理対象ブロック(3×3の大きさ)上で画素データを順次スキャンすることにより、急激な画素データの変化点があるか否かを検出する。この判断の結果として、画像の急激な変化領域が検出された場合は(S2300:YES)、ステップS2350に進み、画像の急激な変化領域が検出されなかった場合は(S2300:NO)、ステップS2400に進む。
【0121】
具体的な検出方法の一例は、次の通りである。ある一定量の過去の観測量の平均値y’(n−1)と現在の観測量y(n)との差の絶対値、つまり、|y’(n−1)−y(n)|の値をしきい値αと比較する。そして、その値がしきい値α以上の場合、つまり、|y’(n−1)−y(n)|≧αの場合は、急激な画素データの変化点があると判断し、その値がしきい値α未満の場合、つまり、|y’(n−1)−y(n)|<αの場合は、急激な画素データの変化点はないと判断する。例えば、図26において、時刻n−1の処理対象ブロック301における所定部分302の観測量の平均値y’(n−1)と、時刻nの処理対象ブロック303における注目画素の部分304の観測量y(n)とを比較する。
【0122】
ステップS2350では、ステップS2250で選択した処理対象ブロックに対して、発明手法による画像復元処理を行う。発明手法による画像復元処理は、上記のように、急激な画素データの変化点がある場合、つまり、エッジ部分などを含む場合であっても、高精度な復元処理が可能である。なお、発明手法による画像復元処理の詳細な手順の一例については、図13のフローチャートを用いて説明した通りである。
【0123】
一方、ステップS2400では、ステップS2250で選択した処理対象ブロックに対して、発明手法以外の他の画像復元手法による画像復元処理を行う。他の画像復元手法としては、ウィナーフィルタを用いた画像復元手法や射影フィルタを用いた画像復元手法、カルマンフィルタを用いた画像復元手法などを含む任意の画像復元手法を用いることができる。急激な画素データの変化点がない場合、つまり、エッジ部分などを含まない場合には、他の画像復元手法でも精度良く復元処理を行うことができるからである。
【0124】
そして、ステップS2450では、ステップS2350での復元処理結果またはステップS2400での復元処理結果を、順次、記憶部150(主記憶装置152または補助記憶装置154)の作業用メモリに保存する。
【0125】
そして、ステップS2500では、時刻nのカウンタの値を1だけインクリメントする。
【0126】
そして、ステップS2550では、1頁分の画像復元処理が終了したか否かを判断する。この判断は、時刻nのカウンタの値に基づいて行われる。この判断の結果として、1頁分の画像復元処理が終了した場合は(S2550:YES)、ステップS2600に進み、1頁分の画像復元処理が終了していない場合は(S2550:NO)、ステップS2250に戻る。
【0127】
ステップS2600では、1頁分の画像復元処理が終了したため、表示メモリ156の更新処理を行う。すなわち、ここでは、1頁分の画像復元処理が終了した時点で表示メモリ156を更新するように構成されている。
【0128】
そして、ステップS2650では、ステップS2000で取り込んだ画像データの全頁分に対して画像復元処理が終了したか否かを判断する。この判断の結果として、全頁分の画像復元処理が終了していない場合は(S2650:NO)、ステップS2250に戻り、全頁分の画像復元処理が終了した場合は(S2650:YES)、上記一連の処理を終了する。
【0129】
次に、ぼけ情報である点拡がり関数(PSF)に対応する観測遷移行列Mp1の設定と実装の具体例について、図27および図28を用いて説明する。
【0130】
図27は、復元モードの実装の一例を説明するための図である。
【0131】
上記のように、観測遷移行列Mp1は、式(4)で定義された9×9次の行列であって、画像の劣化モデルにおけるぼけの点拡がり関数(PSF)に対応しており、観測遷移行列Mp1を構成する各要素hr,s(r,s=−1,0,1)は、既知であって、予めデータ化して適切に設定されている。すなわち、本実施の形態では、画像によってぼけ具合が変わるため、いわばぼけ情報である観測遷移行列Mp1の各要素hr,sの値の組を復元モードとしていくつか予め設定しておき、使用者は、操作部130の復元モード指定部136を介して、その復元モードを任意に指定できるように構成されている。例えば、図27の例では、デジタルカメラなどによく設定されている夜景モードやスポーツモードなどの撮影モード毎に、予め、その撮影モードに合った適切な観測遷移行列Mp1がそれぞれ設定されている。したがって、使用者は、撮影モードを切り替えることによって、最適な復元モードを使用することができる。なお、復元モード(観測遷移行列Mp1の各要素hr,sの値)自体を自動または手動で再調整することも可能である。
【0132】
図28は、復元モードの実装の他の例を説明するための図である。
【0133】
図28の例では、画像によってぼけ具合が変わるため、まず、図28(A)に示すように、ディスプレイ182上で、処理対象の劣化画像を選択した後、図28(B)に示すように、観測遷移行列Mp1(各要素hr,sの値)を変えた復元画像をいくつかディスプレイ182に表示する。このとき、使用者が、例えば、操作部130(タッチパネルなど)を介して、ディスプレイ182に表示された複数の復元画像の中から、最適な復元画像を選択すると、復元画像が決定される。なお、使用者の選択の履歴をデータとして蓄積することによって、使用者の好みに合った復元モード(観測遷移行列Mp1の各要素hr,sの値)を自動で選択することも可能となる。
【0134】
このように、本実施の形態によれば、ARシステムのコンセプトを必要としない新しい状態空間モデル(状態方程式と観測方程式)を構成して1段階処理の新しい予測法によって画像復元を行うため、シンプルで実用的かつ高性能な画像復元を実現することができる。すなわち、AR次数の決定およびAR係数を推定するステップを必要としないシンプルな構成で、かつ、非定常性が強い自然画像に対しても有効に復元処理を行うことができる実用性を有し、しかも、従来手法1〜3に比べて画像復元能力を向上することができる。
【0135】
また、本実施の形態では、上記のように、ARシステムのコンセプトを必要としないため、画像の撮り直しが許されない場面にも広く適用可能である。なお、もちろん、ここでも、上記のように、画像の種類は限定されず、静止画像でも動画像でもよい。
【0136】
例えば、古文書の傷や汚れなどをぼけや雑音とみなし、発明手法1を用いてこれらを取り除き、文字・物体認識技術を用いることによって、古文書の解読、復元、およびデジタルデータベース化に効果を発揮することが期待される。
【0137】
また、防犯分野においても、近年、防犯カメラ(監視カメラ)などの設置が普及しつつあり、犯罪者の多くは防犯カメラなどに撮影されていることが多いため、この画像が犯罪者の発見の手がかりとなることがある。しかし、これらの画像は一般に画質が荒く、動いている対象物の画像はぼけが生じていることが一般的である。また、これらの画像を拡大した場合には、通常、画像はより一層劣化する。そこで、拡大画像を含めてこれらの劣化画像に発明手法1を適用することによって、劣化画像のぼけと雑音を取り除き、クリアな画像を提供することで、犯罪者の早期発見につながることが期待される。なお、このように、拡大画像を含めて劣化画像に発明手法1を適用することによって、劣化画像のぼけと雑音を取り除き、クリアな画像を提供することは、防犯分野に限らず、監視カメラなどで撮影された画像に基づいて、事故の判定や器具の故障診断などを行う場合にも適用可能である。
【0138】
また、近年、カメラ機能付き携帯電話やデジタルカメラなどが急速な勢いで普及し、1人1台の時代となっている。しかし、これらの製品に搭載されている画像復元技術は、顔認識やフィルタなどのぼけや雑音を予防する技術であり、撮り直しが許されない場面を想定していない。したがって、これらの製品に発明手法1を適用することによって、撮り直しが許されない瞬間的な画像復元が可能になる。
【0139】
また、医療の分野において、健康を調査する最も効果的な手段の1つとして、内視鏡やレントゲン、エコー、CT、MRIなどの画像診断がある。例えば、肺や心臓などの画像診断では、肺や心臓などの動作によってぼけと雑音を含んだ情報から健康診断を行っている。これらの画像診断に発明手法を適用してぼけと雑音を取り除くことによって、専門医にクリアな画像を迅速に提供することが可能となり、病気の早期発見の一助になることが期待される。
【0140】
また、近年カーナビゲーションシステムが普及していることから、フロント・リアカメラを搭載した自動車も増えている。このようなリアルタイムで画像復元処理が必要な自動車のフロント・リアカメラなどの場合においても発明手法1は有効である。
【0141】
なお、上記のように、本発明に係る画像復元方法の具体的な手法(発明手法)は、発明手法1に限定されない。以下では、他の発明手法として、5つの発明手法2〜6を順次説明する。ただし、発明手法1〜6相互間の比較を容易にするため、上記した発明手法1についても簡単にまとめて説明しておく。また、共通する部分については、適宜、その説明を省略する。
【0142】
<発明手法1>
まず、各発明手法1〜6に共通の問題設定として、画像復元の理論となる画像劣化モデルについて、再度、簡単に説明しておく。図29は、劣化画像モデルを視覚的に説明するための図である(図3、式(1)参照)。
【0143】
画像の劣化は、図29の上側に示すように、原画像にぼけと雑音が加わることによって生じる。ぼけは、ある画素がその周囲の画素に影響することによって生じるが、雑音は、画素に関係なく発生する。カメラなどで撮影された画像の場合、ぼけは、手ぶれや焦点ずれなどに起因して生じ、雑音は、暗電流や熱雑音などに起因して不可避的に発生する。劣化画像は、一般的に、原画像とぼけの点拡がり関数(PSF)の畳み込み(コンボリューション)に雑音を加えたモデルによって得られる。具体的には、劣化画像g(x,y)は、原画像をf(m,n)、ぼけの点拡がり関数をh(x,y)、雑音をn(x,y)でそれぞれ表すと、図29に示す劣化画像生成式で表される。この劣化画像モデルでは、雑音n(x,y)と原画像f(m,n)は無相関であること、雑音n(x,y)は加法性白色ガウス雑音であること、が条件となっている。劣化画像の復元は、劣化画像からぼけと雑音を取り除くことにほかならない。
【0144】
図30は、発明手法1の状態空間モデルを視覚的に説明するための図、図31は、発明手法1の状態過程を視覚的に説明するための図、図32は、発明手法1の観測過程を視覚的に説明するための図である(以上、図4〜図11、式(2)〜式(4)参照)。
【0145】
発明手法1では、図30にまとめて示すように、新しい状態空間モデル(状態方程式と観測方程式)を構成し、定式化している。この状態空間モデルは、状態過程と観測過程とで構成される。状態過程は、図31に示すように、状態方程式によって記述され、観測過程は、図32に示すように、観測方程式によって記述される。ここで、状態方程式において、ベクトルxp1(n)は状態ベクトル(原画像情報)、行列Φp1は状態遷移行列、ベクトルδp1(n)は駆動源ベクトルである。また、観測方程式において、ベクトルyp1(n)は観測ベクトル(観測画像情報、劣化画像情報)、行列Mp1は観測遷移行列、ベクトルεp1(n)は観測雑音ベクトルである。また、ベクトルxp1(n)、δp1(n)、yp1(n)、εp1(n)の変数としての「n」は、装置の時刻nである。時刻nは、注目画素(または注目画素領域)と複数の周辺画素(または周辺画素領域)とからなる処理対象ブロックの処理の順番(n番目)を表している。状態方程式は、観測対象のシステムを状態空間モデルで記述したものであり、内部状態つまり状態変数(ここでは、状態ベクトルxp1(n)、xp1(n+1))の時間に対する生成過程を表している。また、観測方程式は、何らかの観測装置を通じて観測する過程を記述したものであり、観測結果(ここでは、観測ベクトルyp1(n))が、被観測量つまり入力(ここでは、状態ベクトルxp1(n))に依存して生成される様子を示している。
【0146】
図33は、発明手法1の状態方程式と観測方程式をまとめて示す説明図であり、図34は、発明手法1のアルゴリズムを説明するための図である(図12、図14、式(5)〜式(10)参照)。なお、式の表記法として、ここでは、一般に、アルファベットの小文字「αγ」は、サイズがγ×1のベクトルを表し、アルファベットの大文字「Βγ」は、サイズがγ×γの行列を表している。
【0147】
図33に示す状態方程式および観測方程式、ならびに図34に示すアルゴリズムにおいて、Kは周辺画素領域(厳密には処理対象ブロック)に含まれる画素数であり、m、nをそれぞれ周辺画素領域の縦軸と横軸の長さとした場合、K=m×nで表される。すなわち、発明手法1では、劣化画像の復元処理は、復元対象である注目画素のみを用いるのではなく、その周辺の画素も含めて処理を行う。具体的には、例えば、J×Jの大きさの画像のある注目画素を中心にして、その周辺m×n(ただし、J>m、n)の領域を処理対象ブロックとし、この処理対象ブロックの中の画素情報を用いて、画像の復元処理を行う。処理対象ブロックは、画像復元においてm×nの画素を用いて処理を行う範囲を意味する。図33および図34の右下には、K=9(=3×3)の場合の処理対象ブロックが例示されている。
【0148】
図35は、図34のアルゴリズムを実行する処理手順を示すフローチャートである(図13参照)。
【0149】
ここでは、図34のアルゴリズムの表記に対応して、行列とベクトルのサイズをそれぞれ明記している。具体的には、n→(n+1)の推定誤差の相関値(行列)(S1100)は、サイズがK×Kの行列、重み係数(行列)(S1200)は、サイズがK×Kの行列、n→(n+1)の状態量の最適推定値(ベクトル)(S1300)は、サイズがK×1のベクトル、(n+1)→(n+1)の状態量の最適推定値(ベクトル)(S1400)は、サイズがK×1のベクトル、(n+1)→(n+1)の推定誤差の相関値(行列)(S1600)は、サイズがK×Kの行列となる。
【0150】
<発明手法2>
図36は、発明手法1の課題を視覚的に説明するための図である。
【0151】
発明手法1では、図36に示すように、時刻nの変化は、画素領域(処理対象ブロック)の遷移を用いて表される。時刻n+1では、例えば、過去の1時刻前の9(=3×3)個の状態量(原画像の画素情報)x〜xを用いて、現在の処理対象ブロックの9(=3×3)個の状態量(原画像の画素情報)x〜x12を推定することにより、注目画素領域の状態量xを復元する。このように、発明手法1では、1つの注目画素領域の状態量(要素)を復元するために、処理対象ブロック内の状態量(要素)をすべて推定している。しかし、そもそも、1つの注目画素領域の要素を復元するためには、その注目画素領域の要素のみを求めればよく、他の要素の推定は不要である。そこで、発明手法2では、注目画素領域のみを考慮した状態空間モデルを再構成する。具体的には、発明手法2では、発明手法1の観測方程式を再構成して注目画素領域の要素のみを抽出するようにしている。これにより、発明手法2では、発明手法1に対して大幅な演算量の軽減を可能にして、高速(低演算量)な劣化画像復元を実現している。
【0152】
図37は、発明手法2の状態方程式の構成を視覚的に示す説明図であり、特に、(A)は、発明手法1の状態方程式の具体例を示す図、(B)は、時刻の変化に伴う処理対象ブロックの変化を示す図である。なお、図37に示すように、発明手法2の状態方程式は、発明手法1の状態方程式と全く同じであるため、その説明を省略する。
【0153】
図38は、発明手法2の観測方程式の構成を視覚的に示す説明図であり、図39は、発明手法2の観測方程式の構成を視覚的に示す他の説明図である。
【0154】
発明手法2では、図38に示すように、発明手法1の観測方程式において注目画素領域の要素(例えば、画素情報y(n))のみを抽出して観測方程式を再構成する。これにより、発明手法2では、図39に示すように、発明手法1における観測ベクトルyp1(n)、観測遷移行列Mp1、観測雑音ベクトルεp1(n)が、それぞれ、スカラー化、ベクトル化、スカラー化されて、観測値yp2(n)、観測遷移ベクトルmp2、観測雑音εp2(n)となる。この結果、発明手法2は、発明手法1に対して、サイズが縮小されるため、演算量が大幅に軽減される。なお、演算量の軽減効果については、後で詳述する。
【0155】
図40は、発明手法2の状態方程式と観測方程式をまとめて示す説明図であり、図41は、発明手法2のアルゴリズムを説明するための図である。なお、式の表記法として、ここでは、一般に、アルファベットの小文字「αγ」は、サイズがγ×1のベクトルを表し、アルファベットの大文字「Βγ」は、サイズがγ×γの行列を表している。
【0156】
図40に示す状態方程式および観測方程式、ならびに図41に示すアルゴリズムにおいて、Kは周辺画素領域(厳密には処理対象ブロック)に含まれる画素数であり、m、nをそれぞれ周辺画素領域の縦軸と横軸の長さとした場合、K=m×nで表される。発明手法2のアルゴリズムは、図42に示すように、初期設定(Initialization)の過程と反復(Iteration)の過程とに大別され、反復過程は、発明手法1の場合(図34参照)よりも演算量が軽減されるよう、図40に示す新しい状態空間モデル(状態方程式と観測方程式)を基に構成されている。反復過程では、1〜6の手順を逐次繰り返す。
【0157】
図42は、図41のアルゴリズムを実行する処理手順を示すフローチャートである。
【0158】
上記のように、発明手法2では、図39に示すように、発明手法1における観測ベクトルyp1(n)、観測遷移行列Mp1、観測雑音ベクトルεp1(n)が、それぞれ、スカラー化、ベクトル化、スカラー化されて、観測値yp2(n)、観測遷移ベクトルmp2、観測雑音εp2(n)となっている。これにより、図42のフローチャートにおいて、n→(n+1)の推定誤差の相関値(S3100)は、サイズがK×Kの行列、重み係数(S3200)は、サイズがK×1のベクトル、n→(n+1)の状態量の最適推定値(S3300)は、サイズがK×1のベクトル、(n+1)→(n+1)の状態量の最適推定値(S3400)は、サイズがK×1のベクトル、(n+1)→(n+1)の推定誤差の相関値(S3600)は、サイズがK×Kの行列となる。
【0159】
図43は、発明手法1と発明手法2のアルゴリズムを比較した説明図である。
【0160】
図43に示すように、発明手法1のアルゴリズムでは、逆行列の計算が必要である(反復過程つまり更新式の手順2参照)。一般に、逆行列の計算は複雑であり、演算量が多い。これに対し、発明手法2では、上記したサイズの縮小により、発明手法1における逆行列の計算が、スカラーの逆数の計算となる(反復過程つまり更新式の手順2参照)。スカラーの逆数の計算は簡単であり、演算量は1回で済む。したがって、発明手法2は、発明手法1よりも演算量が大幅に軽減される。
【0161】
図44は、発明手法2の演算量評価を説明するための図である。なお、ここで、従来手法3は、カルマンフィルタを用いた画像復元手法(例えば、非特許文献6参照)である。ただし、非特許文献6に記載された手法は、原画像に雑音が加わった劣化画像にしか対応できないため、本発明者は、非特許文献6に記載された手法をぼけと雑音に対応できるように拡張した上で、以下の各種比較(シミュレーションを含む)を行った。
【0162】
図44において、演算量は、乗算の回数によって表される。図44の上側に、従来手法3、発明手法1、および発明手法2に対する、更新1回当たりの演算回数の比較結果を示す。これに基づき、更新1回当たりの演算量(乗算回数)を求めると、図44の中央左側に示す棒グラフが得られる。この棒グラフから明らかなように、発明手法2は、発明手法1よりも演算量が大幅に減少している。演算量の減少効果で見た場合、その大小関係は、従来手法3<発明手法1<発明手法2、となる。したがって、発明手法2は、大幅な演算量の軽減を実現している。
【0163】
図45は、発明手法2の復元処理速度評価を説明するための図である。
【0164】
ここでは、各手法の処理速度を比較するために、図45に示すスペック表に基づいて、サイズが異なる2種類の画像の処理時間を測定した。その結果、発明手法1および発明手法2は、いずれも、従来手法3と比較して、AR係数推定(導出)の演算量が不要であるため、処理速度が速くなっている。また、発明手法2は、上記のように発明手法1よりも演算量が大幅に少ないため、発明手法1よりも高速な処理が可能となっている。処理速度で見た場合、その大小関係は、図45に示すように、従来手法3<発明手法1<発明手法2、となる。したがって、発明手法2は、高速な画像復元処理が可能である。
【0165】
また、本発明者は、発明手法2の有効性を検証するためにシミュレーションを行った。具体的には、発明手法2の画像復元能力を評価するために、視覚評価と客観評価を行った。視覚評価は、原画像と復元画像とを視覚により比較した評価であり、客観評価は、数値による評価である。
【0166】
図46は、シミュレーションの条件を説明するための図である。
【0167】
本シミュレーションでも、図15のシミュレーション条件と同様に、原画像として、「カメラマン(Cameraman)」と「レナ(Lenna)」を用いた(画像サイズは共に256×256)。また、原画像に加えるぼけに対応する点拡がり関数(PSF)のモデルとして、図15に示す2次元ガウス関数を用い、雑音として、加法性白色ガウス雑音を用いた。また、信号対雑音比SNRを20dBとした。
【0168】
図47は、比較対象手法との対比において原画像「カメラマン」に対する発明手法2のシミュレーション結果(視覚評価)を示す図であり、図48は、図47中の破線の円の周辺領域を拡大した図である。また、図49は、比較対象手法との対比において原画像「レナ」に対する発明手法2のシミュレーション結果(視覚評価)を示す図であり、図50は、図49中の破線の円の周辺領域を拡大した図である。なお、ここで、比較対象手法は、従来手法3と発明手法1である。
【0169】
図47〜図50に示すように、発明手法1および発明手法2は、従来手法3に比べて、明らかに原画像を忠実に復元していることがわかる。ただし、発明手法1と発明手法2を比較した場合、画像復元能力はほぼ同等のレベルであるといえる。
【0170】
図51は、比較対象手法との対比において原画像に対する発明手法2のシミュレーション結果(客観評価)を示す図である。
【0171】
ここでは、従来手法3、発明手法1、および発明手法2の画像復元能力を数値により評価するため、図51に示すPSNR(Peak Signal-to-Noise Ratio)(式(11)参照)を用いて画像復元能力を評価した。なお、PSNRは、画像の似かより度を示す値であり、値が大きいほど劣化(ぼけや雑音など)が少なく画像として良好であるといえる。
【0172】
この場合、図51に示すように、「カメラマン」、「レナ」どちらの画像についても、発明手法1および発明手法2は、従来手法3よりもPSNRの数値が大きいことが確認できる。これにより、発明手法1および発明手法2は、客観評価の点でも、従来手法3に比べて画像復元能力が高いことがわかる。一方、発明手法1と発明手法2を比較した場合には、発明手法1よりも発明手法2のほうが若干高い客観評価(数値による評価)となっている。これにより、客観評価の大小(良否)関係は、図51に示すように、従来手法3<発明手法1<発明手法2、となる。ただし、この若干の差異は、図47〜図50の比較からわかるように、ほぼないに等しい。
【0173】
以上より、発明手法2は、発明手法1に対して、画像復元能力をほぼ維持しながら、演算量を大幅に軽減して高速な画像復元処理を実現していることがわかる。
【0174】
<発明手法3>
図52は、発明手法2の課題を視覚的に説明するための図である。
【0175】
発明手法2は、上記のように、発明手法1に対して、演算量を大幅に軽減して高速な画像復元処理を実現することができるものの、画像復元能力が若干低下している。これは、図52に示すように、発明手法2の観測方程式において、有色性駆動源の影響が大きい領域が含まれていないためであると考えられる。そこで、発明手法3では、注目画素領域のみでなく有色駆動源をも考慮した状態空間モデルを再構成する。すなわち、発明手法3では、有色性駆動源の影響を最大限に与えるように状態空間モデルを再構成する。
【0176】
図53は、発明手法3の観測方程式の構成を視覚的に示す説明図であり、図54は、発明手法3の状態空間モデルを視覚的に説明するための図である。なお、発明手法3の状態方程式は、発明手法2(ひいては発明手法1)の状態方程式と全く同じであるため、その説明を省略する。
【0177】
発明手法3では、図53に示すように、発明手法1の観測方程式において、注目画素領域の要素(例えば、画素情報y(n))に加えて、有色性駆動源の影響を受ける領域の要素(例えば、画素情報y(n)、y(n)、y(n))をも抽出して観測方程式を再構成する。すなわち、発明手法3では、有色性駆動源の影響を考慮するために、発明手法2と比較して、画素情報y(n)、y(n)、y(n)をも取り入れる。これにより、発明手法3では、発明手法1における観測ベクトルyp1(n)、観測遷移行列Mp1、および観測雑音ベクトルεp1(n)は、サイズがそれぞれ縮小されて、観測ベクトルyp3(n)、観測遷移行列Mp3、観測雑音ベクトルεp3(n)となる。このように、発明手法3では、図54に示すように、状態方程式における有色性駆動源の影響を最大限に考慮するべく、注目画素領域と有色性駆動源の影響を受ける領域とから観測方程式を再構成する。これにより、発明手法3は、発明手法2の場合よりもサイズの縮小は小さいものの、発明手法1に対してサイズが縮小されるため、演算量は大幅に軽減される。
【0178】
図55は、発明手法3の状態方程式と観測方程式をまとめて示す説明図であり、図56は、発明手法3のアルゴリズムを説明するための図である。なお、式の表記法として、ここでは、一般に、アルファベットの小文字「αγ」は、サイズがγ×1のベクトルを表し、アルファベットの大文字「Βγ」は、サイズがγ×γの行列を表し、アルファベットの大文字「C(γ×ζ)」は、サイズがγ×ζの行列を表している。
【0179】
図55に示す状態方程式および観測方程式、ならびに図56に示すアルゴリズムにおいて、Kは周辺画素領域(厳密には処理対象ブロック)に含まれる画素数であり、m、nをそれぞれ周辺画素領域の縦軸と横軸の長さとした場合、K=m×nで表される。また、Lは、注目画素領域と有色性駆動源を含む領域とを合わせた領域に含まれる画素数であり、L=1+mで表される(L<K)。発明手法3のアルゴリズムは、図56に示すように、初期設定(Initialization)の過程と反復(Iteration)の過程とに大別され、反復過程は、発明手法2の場合(図41参照)よりも演算量が若干増加するものの、発明手法1の場合(図34参照)よりも演算量が軽減されるよう、図55に示す新しい状態空間モデル(状態方程式と観測方程式)を基に構成されている。反復過程では、1〜6の手順を逐次繰り返す。
【0180】
図57は、図56のアルゴリズムを実行する処理手順を示すフローチャートである。
【0181】
上記のように、発明手法3では、図53に示すように、発明手法1における観測ベクトルyp1(n)、観測遷移行列Mp1、観測雑音ベクトルεp1(n)は、サイズがそれぞれ縮小されて(つまり、L<K)、観測ベクトルyp3(n)、観測遷移行列Mp3、観測雑音ベクトルεp3(n)となっている。これにより、図57のフローチャートにおいて、n→(n+1)の推定誤差の相関値(S4100)は、サイズがK×Kの行列、重み係数(S4200)は、サイズがK×Lの行列、n→(n+1)の状態量の最適推定値(S4300)は、サイズがK×1のベクトル、(n+1)→(n+1)の状態量の最適推定値(S4400)は、サイズがK×1のベクトル、(n+1)→(n+1)の推定誤差の相関値(S4600)は、サイズがK×Kの行列となる。
【0182】
図58は、発明手法3の演算量評価を説明するための図である。
【0183】
図58の上側に、発明手法3に対する、更新1回当たりの演算回数(乗算回数)を示す。これに基づき、更新1回当たりの演算量(乗算回数)を求め、得られた結果を図44に示す棒グラフに加えると、図58の中央左側に示す棒グラフが得られる。この棒グラフから明らかなように、発明手法3は、発明手法1よりも演算量は大幅に少ない。演算量の減少効果で見た場合、その大小関係は、従来手法3<発明手法1<発明手法3<発明手法2、となる。したがって、発明手法3は、良好な演算量軽減を実現している。
【0184】
図59は、発明手法3の復元処理速度評価を説明するための図である。
【0185】
ここでは、図45の結果に、発明手法3の測定結果を追加した。その結果、発明手法3も、従来手法3と比較して、AR係数推定(導出)の演算量が不要であるため、処理速度が速くなっている。また、発明手法3は、上記のように発明手法1よりも演算量が大幅に少ないため、発明手法1よりも高速な処理が可能となっている。発明手法3は、発明手法2よりも処理速度が若干遅い程度である。処理速度で見た場合、その大小関係は、図59に示すように、従来手法3<発明手法1<発明手法3<発明手法2、となる。したがって、発明手法3も、高速な画像復元処理が可能である。
【0186】
図60は、比較対象手法との対比において原画像「カメラマン」に対する発明手法3のシミュレーション結果(視覚評価)を示す図であり、図61は、図60中の破線の円の周辺領域を拡大した図である。また、図62は、比較対象手法との対比において原画像「レナ」に対する発明手法3のシミュレーション結果(視覚評価)を示す図であり、図63は、図62中の破線の円の周辺領域を拡大した図である。なお、シミュレーション条件は、図46に示す通りである。また、ここでの比較対象手法は、従来手法3と発明手法1である。
【0187】
図60〜図63に示すように、発明手法3は、従来手法3に比べて、明らかに原画像を忠実に復元していることがわかる。しかも、発明手法1と比較しても、画像復元能力は、発明手法3のほうが高いことがわかる。
【0188】
図64は、比較対象手法との対比において原画像に対する発明手法3のシミュレーション結果(客観評価)を示す図である。
【0189】
この場合、図64に示すように、「カメラマン」、「レナ」どちらの画像についても、発明手法3は、他のどの手法(従来手法3はもとより発明手法1および発明手法2)よりもPSNRの数値が大きいことが確認できる。これにより、発明手法3は、客観評価の点でも、他の手法(従来手法3、発明手法1、発明手法2)に比べて画像復元能力が高いことがわかる。これにより、客観評価(数値による評価)の大小(良否)関係は、図64に示すように、従来手法3<発明手法1<発明手法2<発明手法3、となる。
【0190】
以上より、発明手法3は、発明手法2よりも演算量を若干増加させつつも、発明手法1および発明手法2よりも高い良好な画像復元能力を発揮することができる高速な画像復元処理を実現していることがわかる。
【0191】
<発明手法4>
図65は、1次元の画像復元手法の特徴を説明するための図である。
【0192】
発明手法1〜3は、1次元の画像復元手法である。1次元の画像復元手法である発明手法1〜3では、図65に示すように、周辺画素領域を用いて注目画素領域(1画素のみ)を復元する。このとき、画像復元手法は、用途や状況などに応じて適用可能である。例えば、上記のように、発明手法2は、発明手法1に比べて、高速(低演算量)な画像復元手法であり、発明手法3は、発明手法2よりも演算量は若干多いが、発明手法1および発明手法2よりも画像復元能力を向上させた手法であるため、処理速度を優先したい場合には発明手法2を適用し、画質を優先したい場合には発明手法3を適用する、という使い方が可能である。
【0193】
図66は、2次元の画像復元手法の特徴を説明するための図である。
【0194】
発明手法4では、1次元の発明手法1を2次元に拡張した2次元の画像復元手法を提案する。発明手法4では、発明手法1の場合よりも注目画素領域を大きくとることにより、更新回数を軽減して処理効率のよい復元が可能となる。そのため、発明手法4では、発明手法1における画素領域(1画素のみ)を、縦と横の任意の長さをp、qとした画素領域に拡張する(以下「画素領域の2次元化」ともいう)。1つの画素領域がp=q=2、かつ、処理対象ブロック(注目画素領域と周辺画素領域とからなる)が3pq×3pqの場合の例を、図66に示す。画素領域の2次元化に当たり、画素領域のサイズおよび処理対象ブロックのサイズは、それぞれ、自由に選択可能である。画素領域の2次元化により、注目画素領域の画素数が多くなるため、処理の効率化が図られる。また、p=q=1の場合、注目画素領域は1次元の画像復元手法(発明手法1〜3)と等しくなるため、2次元の画像復元手法である発明手法4は、柔軟性のある手法といえる。なお、上記のように、1次元カルマンフィルタを2次元に拡張した2次元カルマンフィルタも知られているが、2次元カルマンフィルタも、1次元カルマンフィルタと同様に、AR係数を用いるため、AR係数の推定精度に起因する画像復元能力の低下は避けられない。これを解決するためには、AR係数を必要としない新しい予測法を発案して、2次元カルマンフィルタの状態空間モデルを再構成すればよい。発明手法4は、まさに、このための新しい手法であり、発明手法1を2次元に拡張した新しい予測法に基づく画像復元手法である。以下、発明手法4について、より詳細に説明する。
【0195】
図67は、発明手法4の状態空間モデルを視覚的に説明するための図である。また、図68は、発明手法4の状態過程を視覚的に説明するための図であり、図69は、発明手法4の状態方程式の構成を視覚的に示す説明図である。
【0196】
状態過程について、発明手法4では、図68に示すように、画素領域として、サイズがp×qの領域を指定し、指定した画素領域に含まれる画素情報を並べてpq×1次の画素領域ベクトルを生成する。例えば、図68に示すように、p=q=2、かつ、3pq×3pqの領域の場合を例にとると、一番左上の画素領域ベクトルx(n)は、4個の画素情報x(n)、x(n)、x(n)、x(n)を要素とする4×1次のベクトルであり、一番右下の画素領域ベクトルx(n)は、4個の画素情報x33(n)、x34(n)、x35(n)、x36(n)を要素とする4×1次のベクトルである。そして、9個の画素領域に対応する画素領域ベクトルx(n)(k=1、2、…、9)を並べて状態ベクトルx(n)を表現する。すなわち、発明手法4では、図69に示すように、画素領域を、発明手法1における画素情報と対応付けて、同様の処理を行う。
【0197】
図70は、発明手法4の観測過程を視覚的に説明するための図であり、図71は、発明手法4の観測方程式の構成を視覚的に示す説明図である。
【0198】
観測過程についても、発明手法4では、図70に示すように、画素領域として、サイズがp×qの領域を指定し、指定した画素領域に含まれる画素情報を並べてpq×1次の画素領域ベクトルを生成する。例えば、図70に示すように、p=q=2、かつ、3pq×3pqの領域の場合を例にとると、一番左上の画素領域ベクトルy(n)は、4個の画素情報y(n)、y(n)、y(n)、y(n)を要素とする4×1次のベクトルであり、一番右下の画素領域ベクトルy(n)は、4個の画素情報y33(n)、y34(n)、y35(n)、y36(n)を要素とする4×1次のベクトルである。そして、9個の画素領域に対応する画素領域ベクトルy(n)(k=1、2、…、9)を並べて観測ベクトルy(n)を表現する。すなわち、発明手法4では、図71に示すように、画素領域を、発明手法1における画素情報と対応付けて、同様の処理を行う。
【0199】
図72は、発明手法4の観測遷移行列の、p=q=2、かつ、3pq×3pqの領域の場合の具体例を示す図であり、図73は、図72の観測遷移行列の構成を視覚的に示す説明図である。また、図74は、図72の観測遷移行列を構成する係数の割り当て方法の具体例を示す図、図75は、図72の観測遷移行列を構成する係数の割り当て方法の具体例を示す他の図、図76は、図72の観測遷移行列を構成する係数の割り当て方法の具体例を示すさらに他の図である。
【0200】
ここでは、発明手法4の観測遷移行列を構成する3×3(=9)個の係数hr,s(座標r,s=−1,0,1)は、発明手法1の場合と同じであり、画素情報を単位として見たときに、発明手法1の場合と同様の関係が成り立つように、割り当てられている。係数hr,s(r,s=−1,0,1)の具体的な割り当て方法は、画素情報を単位としたときの観測量yの位置に応じて、図74〜図76にそれぞれ示す通りである。
【0201】
図77は、発明手法4の状態方程式と観測方程式をまとめて示す説明図であり、図78は、発明手法4のアルゴリズムを説明するための図である。なお、式の表記法として、ここでは、一般に、アルファベットの小文字「αγ」は、サイズがγ×1のベクトルを表し、アルファベットの大文字「Βγ」は、サイズがγ×γの行列を表している。
【0202】
図77に示す状態方程式および観測方程式、ならびに図78に示すアルゴリズムにおいて、Kは周辺画素領域(厳密には処理対象ブロック)に含まれる画素数であり、m、nをそれぞれ周辺画素領域の縦軸と横軸の長さとした場合、K=m×nで表される。ここで、m、nは、それぞれ、m=l×p、n=l×qで表される。p、qは、それぞれ、注目画素領域の縦軸と横軸の長さであり、l、lは、それぞれ、p×qと同等のサイズを周辺画素領域の縦軸と横軸に構成した数である。
【0203】
図79は、図78のアルゴリズムを実行する処理手順を示すフローチャートである。
【0204】
ここでは、図78のアルゴリズムの表記に対応して、行列とベクトルのサイズをそれぞれ明記している。具体的には、n→(n+1)の推定誤差の相関値(行列)(S5100)は、サイズがK×Kの行列、重み係数(行列)(S5200)は、サイズがK×Kの行列、n→(n+1)の状態量の最適推定値(ベクトル)(S5300)は、サイズがK×1のベクトル、(n+1)→(n+1)の状態量の最適推定値(ベクトル)(S5400)は、サイズがK×1のベクトル、(n+1)→(n+1)の推定誤差の相関値(行列)(S5600)は、サイズがK×Kの行列となる。
【0205】
図80は、発明手法4の演算量評価を説明するための図である。なお、ここで、従来手法4は、2次元カルマンフィルタを用いた画像復元手法(例えば、非特許文献8参照)である。ただし、非特許文献8に記載された手法は、原画像に雑音が加わった劣化画像にしか対応できないため、本発明者は、非特許文献8に記載された手法をぼけと雑音に対応できるように拡張した上で、以下の各種比較(シミュレーションを含む)を行った。
【0206】
図80の上側に、従来手法4および発明手法4に対する、更新1回当たりの演算回数(乗算回数)を示す。これに基づき、256×256サイズの1枚の画像を復元するのに必要な演算量(乗算回数)を求め、得られた結果を発明手法1に対する結果と一緒に示すと、図80の中央左側に示す棒グラフが得られる。この棒グラフから明らかなように、発明手法4は、従来手法4よりも演算量は少ないが、発明手法1よりも演算量は多い。しかし、発明手法4は、発明手法1に比べて、注目画素領域が広いため、更新回数は少ない。演算量の減少効果で見た場合、その大小関係は、従来手法4<発明手法4<発明手法1、となる。したがって、発明手法4は、注目画素領域と周辺画素領域のサイズをそれぞれ変化させることにより、効率のよい画像復元が可能である。
【0207】
図81は、発明手法4の復元処理速度評価を説明するための図である。
【0208】
ここでは、発明手法1の測定結果に、従来手法4と発明手法4の測定結果を追加した。その結果、発明手法4は、従来手法4と比較して、AR係数推定(導出)の演算量が不要であるため、処理速度が速くなっている。また、発明手法4は、上記のように発明手法1よりも演算量は多いものの更新回数は少ないため、発明手法1よりも高速な処理が可能となっている。処理速度で見た場合、その大小関係は、図81に示すように、従来手法4<発明手法1<発明手法4、となる。したがって、発明手法4は、高速な画像復元処理が可能である。
【0209】
図82は、比較対象手法との対比において原画像「カメラマン」に対する発明手法4のシミュレーション結果(視覚評価)を示す図であり、図83は、図82中の破線の円の周辺領域を拡大した図である。また、図84は、比較対象手法との対比において原画像「レナ」に対する発明手法4のシミュレーション結果(視覚評価)を示す図であり、図85は、図84中の破線の円の周辺領域を拡大した図である。なお、シミュレーション条件は、図46に示す通りである。また、ここでの比較対象手法は、従来手法3と発明手法1である。
【0210】
図82〜図85に示すように、発明手法4は、従来手法3に比べて、明らかに原画像を忠実に復元していることがわかる。しかも、発明手法1と比較しても、画像復元能力は、発明手法4のほうが高いことがわかる。
【0211】
図86は、比較対象手法との対比において原画像に対する発明手法4のシミュレーション結果(客観評価)を示す図である。
【0212】
この場合、図86に示すように、「カメラマン」、「レナ」どちらの画像についても、発明手法4は、発明手法3よりもPSNRの数値はやや小さいものの、従来手法3はもとより発明手法1および発明手法2よりもPSNRの数値が大きいことが確認できる。これにより、発明手法4は、客観評価の点でも、従来手法3、発明手法1、および発明手法2に比べて画像復元能力が高いことがわかる。これにより、客観評価(数値による評価)の大小(良否)関係は、図86に示すように、従来手法3<発明手法1<発明手法2<発明手法4<発明手法3、となる。
【0213】
以上より、発明手法4は、発明手法1よりも演算量は多いものの、発明手法1および発明手法2よりも高い良好な画像復元能力を発揮することができる高速な画像復元処理を実現していることがわかる。
【0214】
<発明手法5>
図87は、発明手法4の課題を視覚的に説明するための図である。
【0215】
発明手法4では、図87に示すように、時刻nの変化は、画素領域(処理対象ブロック)の遷移を用いて表される。時刻n+1では、例えば、過去の1時刻前の9(=3×3)個の画素領域にそれぞれ含まれる各状態量(原画像の画素情報)を用いて、現在の処理対象ブロックの9(=3×3)個の画素領域にそれぞれ含まれる各状態量(原画像の画素情報)を推定することにより、注目画素領域の状態量x29、x30、x31、x32を復元する。このように、発明手法4では、1つの注目画素領域に含まれる各状態量(要素)を復元するために、処理対象ブロック内の状態量(要素)をすべて推定している。しかし、そもそも、1つの注目画素領域内の要素を復元するためには、その注目画素領域内の要素のみを求めればよく、他の要素の推定は不要である。そこで、発明手法5では、注目画素領域のみを考慮した状態空間モデルを再構成する。具体的には、発明手法5では、発明手法4の観測方程式を再構成して注目画素領域に含まれる要素のみを抽出するようにしている。これにより、発明手法5では、発明手法4に対して大幅な演算量の軽減を可能にして、高速(低演算量)な劣化画像復元を実現している。なお、発明手法5は、発明手法2を2次元に拡張したものともいえる。
【0216】
図88は、発明手法5の状態方程式の構成を視覚的に示す説明図であり、特に、(A)は、発明手法4の状態方程式の具体例を示す図、(B)は、時刻の変化に伴う処理対象ブロックの変化を示す図である。なお、図88に示すように、発明手法5の状態方程式は、発明手法4の状態方程式と全く同じであるため、その説明を省略する。
【0217】
図89は、発明手法5の観測方程式の構成を視覚的に示す説明図であり、図90は、発明手法5の観測方程式の構成を視覚的に示す他の説明図である。
【0218】
発明手法5では、図89に示すように、発明手法4の観測方程式において注目画素領域に含まれる要素(例えば、画素領域ベクトルy(n))のみを抽出して観測方程式を再構成する。これにより、発明手法5では、図90に示すように、発明手法4における観測ベクトルyp4(n)、観測遷移行列Mp4、観測雑音ベクトルεp4(n)は、それぞれ、サイズが縮小されて、観測ベクトルyp5(n)、観測遷移行列Mp5、観測雑音ベクトルεp5(n)となる。この結果、発明手法5は、発明手法4に対して、サイズが縮小されるため、演算量が大幅に軽減される。なお、演算量の軽減効果については、後で詳述する。
【0219】
図91は、発明手法5の状態方程式と観測方程式をまとめて示す説明図であり、図92は、発明手法5のアルゴリズムを説明するための図である。なお、式の表記法として、ここでは、一般に、アルファベットの小文字「αγ」は、サイズがγ×1のベクトルを表し、アルファベットの大文字「Βγ」は、サイズがγ×γの行列を表し、アルファベットの大文字「C(γ×ζ)」は、サイズがγ×ζの行列を表している。
【0220】
図91に示す状態方程式および観測方程式、ならびに図92に示すアルゴリズムにおいて、Kは周辺画素領域(厳密には処理対象ブロック)に含まれる画素数であり、m、nをそれぞれ周辺画素領域の縦軸と横軸の長さとした場合、K=m×nで表される。ここで、m、nは、それぞれ、m=l×p、n=l×qで表される。p、qは、それぞれ、注目画素領域の縦軸と横軸の長さであり、l、lは、それぞれ、p×qと同等のサイズを周辺画素領域の縦軸と横軸に構成した数である。また、Lは、注目画素領域に含まれる画素数であり、L=p×qで表される。
【0221】
図93は、図92のアルゴリズムを実行する処理手順を示すフローチャートである。
【0222】
上記のように、発明手法5では、図90に示すように、発明手法4における観測ベクトルyp4(n)、観測遷移行列Mp4、観測雑音ベクトルεp4(n)は、それぞれ、サイズが9pq×1→pq×1、9pq×9pq→pq×9pq、9pq×1→pq×1に縮小されて(つまり、L<K)、観測ベクトルyp5(n)、観測遷移行列Mp5、観測雑音ベクトルεp5(n)となっている。これにより、図93のフローチャートにおいて、n→(n+1)の推定誤差の相関値(S6100)は、サイズがK×Kの行列、重み係数(S6200)は、サイズがK×Lの行列、n→(n+1)の状態量の最適推定値(S6300)は、サイズがK×1のベクトル、(n+1)→(n+1)の状態量の最適推定値(S6400)は、サイズがK×1のベクトル、(n+1)→(n+1)の推定誤差の相関値(S6600)は、サイズがK×Kの行列となる。
【0223】
図94は、発明手法4と発明手法5のアルゴリズムを比較した説明図である。
【0224】
図94に示すように、発明手法4のアルゴリズムでは、サイズが9pq×9pqの逆行列の計算が必要であり(反復過程つまり更新式の手順2参照)、演算量が多い。これに対し、発明手法5では、上記したサイズの縮小により、サイズが9pq×9pqの逆行列が、サイズがpq×pqの逆行列となっている(反復過程つまり更新式の手順2参照)。したがって、発明手法5は、発明手法4よりも演算量が大幅に軽減される。
【0225】
図95は、発明手法5の演算量評価を説明するための図である。
【0226】
図95の上側に、発明手法5に対する、更新1回当たりの演算回数(乗算回数)を示す。これに基づき、256×256サイズの1枚の画像を復元するのに必要な演算量(乗算回数)を求め、得られた結果を図80に示す棒グラフに加えると、図95の中央左側に示す棒グラフが得られる。この棒グラフから明らかなように、発明手法5は、従来手法4および発明手法4よりも演算量は大幅に少ない。演算量の減少効果で見た場合、その大小関係は、従来手法4<発明手法4<発明手法1<発明手法5、となる。したがって、発明手法5は、大幅な演算量の軽減が可能であり、高速(低演算量)な画像復元が可能である。
【0227】
図96は、発明手法5の復元処理速度評価を説明するための図
【0228】
ここでは、図81の結果に、発明手法5の測定結果を追加した。その結果、発明手法5も、従来手法4と比較して、AR係数推定(導出)の演算量が不要であるため、処理速度が速くなっている。また、発明手法5は、上記のように発明手法1および発明手法4よりも演算量が少ないため、発明手法1および発明手法4よりも高速な処理が可能となっている。処理速度で見た場合、その大小関係は、図96に示すように、従来手法4<発明手法1<発明手法4<発明手法5、となる。したがって、発明手法5は、高速な画像復元処理が可能である。
【0229】
図97は、比較対象手法との対比において原画像「カメラマン」に対する発明手法5のシミュレーション結果(視覚評価)を示す図であり、図98は、図97中の破線の円の周辺領域を拡大した図である。また、図99は、比較対象手法との対比において原画像「レナ」に対する発明手法5のシミュレーション結果(視覚評価)を示す図であり、図100は、図99中の破線の円の周辺領域を拡大した図である。なお、シミュレーション条件は、図46に示す通りである。また、ここでの比較対象手法は、従来手法3と発明手法1である。
【0230】
図97〜図100に示すように、発明手法5は、従来手法3に比べて、明らかに原画像を忠実に復元していることがわかる。ただし、発明手法1と発明手法5を比較した場合、画像復元能力はほぼ同等のレベルである。
【0231】
図101は、比較対象手法との対比において原画像に対する発明手法5のシミュレーション結果(客観評価)を示す図である。
【0232】
この場合、図101に示すように、「カメラマン」、「レナ」どちらの画像についても、発明手法5は、発明手法3よりもPSNRの数値が小さく、従来手法3よりはPSNRの数値が大きいものの、発明手法1、2、4とはPSNRの数値がほぼ同レベルであることが確認できる。これにより、発明手法5は、客観評価の点でも、従来手法3に比べて画像復元能力が高いことがわかる。これにより、客観評価(数値による評価)の大小(良否)関係は、図101に示すように、従来手法3<発明手法1<発明手法5<発明手法2<発明手法4<発明手法3、となる。
【0233】
以上より、発明手法5は、2次元に拡張した発明手法4とほぼ同等な画像復元能力を維持しつつ、演算量の軽減により高速で、かつ、処理効率が向上した画像復元手法を実現していることがわかる。
【0234】
<発明手法6>
図102は、発明手法5の課題を視覚的に説明するための図である。
【0235】
発明手法5は、上記のように、発明手法4に対して、演算量を大幅に軽減して高速な画像復元処理を実現することができるものの、画像復元能力が若干低下している。これは、図102に示すように、発明手法5の観測方程式において、有色性駆動源の影響が大きい領域が含まれていないためであると考えられる。そこで、発明手法6では、注目画素領域のみでなく有色駆動源をも考慮した状態空間モデルを再構成する。すなわち、発明手法6では、有色性駆動源の影響を最大限に与えるように状態空間モデルを再構成する。
【0236】
図103は、発明手法6の観測方程式の構成を視覚的に示す説明図であり、図104は、発明手法6の状態空間モデルを視覚的に説明するための図である。なお、発明手法6の状態方程式は、発明手法5(ひいては発明手法4)の状態方程式と全く同じであるため、その説明を省略する。
【0237】
発明手法6では、図103に示すように、発明手法4の観測方程式において、注目画素領域に含まれる要素(例えば、画素領域ベクトルy(n))に加えて、有色性駆動源の影響を受ける領域に含まれる要素(例えば、画素領域ベクトルy(n)、y(n)、y(n))をも抽出して観測方程式を再構成する。すなわち、発明手法6では、有色性駆動源の影響を考慮するために、発明手法5と比較して、画素領域ベクトルy(n)、y(n)、y(n)をも取り入れる。これにより、発明手法6では、発明手法4における観測ベクトルyp4(n)、観測遷移行列Mp4、および観測雑音ベクトルεp4(n)は、サイズがそれぞれ縮小されて、観測ベクトルyp6(n)、観測遷移行列Mp6、観測雑音ベクトルεp6(n)となる。このように、発明手法6では、図104に示すように、状態方程式における有色性駆動源の影響を最大限に考慮するべく、注目画素領域と有色性駆動源の影響を受ける領域とから観測方程式を再構成する。これにより、発明手法6は、発明手法5の場合よりもサイズの縮小は小さいものの、発明手法4に対してサイズが縮小されるため、演算量は大幅に軽減される。
【0238】
図105は、発明手法6の状態方程式と観測方程式をまとめて示す説明図であり、図106は、発明手法6のアルゴリズムを説明するための図である。なお、式の表記法として、ここでは、一般に、アルファベットの小文字「αγ」は、サイズがγ×1のベクトルを表し、アルファベットの大文字「Βγ」は、サイズがγ×γの行列を表し、アルファベットの大文字「C(γ×ζ)」は、サイズがγ×ζの行列を表している。
【0239】
図105に示す状態方程式および観測方程式、ならびに図106に示すアルゴリズムにおいて、Kは周辺画素領域(厳密には処理対象ブロック)に含まれる画素数であり、m、nをそれぞれ周辺画素領域の縦軸と横軸の長さとした場合、K=m×nで表される。ここで、m、nは、それぞれ、m=l×p、n=l×qで表される。p、qは、それぞれ、注目画素領域の縦軸と横軸の長さであり、l、lは、それぞれ、p×qと同等のサイズを周辺画素領域の縦軸と横軸に構成した数である。また、Lは、注目画素領域と有色性駆動源を含む領域とを合わせた領域に含まれる画素数であり、L=p×q+m×qで表される(L<K)。発明手法6のアルゴリズムは、図106に示すように、初期設定(Initialization)の過程と反復(Iteration)の過程とに大別され、反復過程は、発明手法5の場合(図92参照)よりも演算量が若干増加するものの、発明手法4の場合(図78参照)よりも演算量が軽減されるよう、図105に示す新しい状態空間モデル(状態方程式と観測方程式)を基に構成されている。反復過程では、1〜6の手順を逐次繰り返す。
【0240】
図107は、図106のアルゴリズムを実行する処理手順を示すフローチャートである。
【0241】
上記のように、発明手法6では、図103に示すように、発明手法4における観測ベクトルyp4(n)、観測遷移行列Mp4、観測雑音ベクトルεp4(n)は、サイズがそれぞれ縮小されて(つまり、L<K)、観測ベクトルyp6(n)、観測遷移行列Mp6、観測雑音ベクトルεp6(n)となっている。これにより、図107のフローチャートにおいて、n→(n+1)の推定誤差の相関値(S7100)は、サイズがK×Kの行列、重み係数(S7200)は、サイズがK×Lの行列、n→(n+1)の状態量の最適推定値(S7300)は、サイズがK×1のベクトル、(n+1)→(n+1)の状態量の最適推定値(S7400)は、サイズがK×1のベクトル、(n+1)→(n+1)の推定誤差の相関値(S7600)は、サイズがK×Kの行列となる。
【0242】
図108は、発明手法6の演算量評価を説明するための図である。
【0243】
図108の上側に、発明手法6に対する、256×256サイズの1枚の画像を復元するのに必要な演算量(乗算回数)示す。これに基づき、更新1回当たりの演算量(乗算回数)を求め、得られた結果を図95に示す棒グラフに加えると、図108の中央左側に示す棒グラフが得られる。この棒グラフから明らかなように、発明手法6は、発明手法5よりも演算量は若干多いものの、発明手法4よりも演算量は大幅に少ない。演算量の減少効果で見た場合、その大小関係は、従来手法4<発明手法4<発明手法1<発明手法6<発明手法5、となる。したがって、発明手法6は、良好な演算量軽減を実現している。
【0244】
図109は、発明手法6の復元処理速度評価を説明するための図である。
【0245】
ここでは、図96の結果に、発明手法6の測定結果を追加した。その結果、発明手法6も、従来手法4と比較して、AR係数推定(導出)の演算量が不要であるため、処理速度が速くなっている。また、発明手法6は、上記のように発明手法4よりも演算量が大幅に少ないため、発明手法4よりも高速な処理が可能となっている。発明手法6は、発明手法5よりも処理速度が若干遅い程度である。処理速度で見た場合、その大小関係は、図109に示すように、従来手法4<発明手法1<発明手法4<発明手法6<発明手法5、となる。したがって、発明手法6も、高速な画像復元処理が可能である。
【0246】
図110は、比較対象手法との対比において原画像「カメラマン」に対する発明手法6のシミュレーション結果(視覚評価)を示す図であり、図111は、図110中の破線の円の周辺領域を拡大した図である。また、図112は、比較対象手法との対比において原画像「レナ」に対する発明手法6のシミュレーション結果(視覚評価)を示す図であり、図113は、図112中の破線の円の周辺領域を拡大した図である。なお、シミュレーション条件は、図46に示す通りである。また、ここでの比較対象手法は、従来手法3と発明手法1である。
【0247】
図110〜図113に示すように、発明手法6は、従来手法3に比べて、明らかに原画像を忠実に復元していることがわかる。しかも、発明手法1と比較しても、画像復元能力は、発明手法6のほうが高いことがわかる。
【0248】
図114は、比較対象手法との対比において原画像に対する発明手法6のシミュレーション結果(客観評価)を示す図である。
【0249】
この場合、図114に示すように、「カメラマン」、「レナ」どちらの画像についても、発明手法6は、他のどの手法(従来手法3はもとより発明手法1〜5)よりもPSNRの数値が大きいことが確認できる。これにより、発明手法6は、客観評価の点でも、他の手法(従来手法3、発明手法1〜5)に比べて画像復元能力が高いことがわかる。これにより、客観評価(数値による評価)の大小(良否)関係は、図114に示すように、従来手法3<発明手法1<発明手法5<発明手法2<発明手法4<発明手法3<発明手法6、となる。
【0250】
以上より、発明手法6は、発明手法5よりも演算量を若干増加させつつも、発明手法1〜5よりも高い良好な画像復元能力を発揮することができる高速な画像復元処理を実現していることがわかる。
【0251】
ここで、発明手法1〜6の状態空間モデル(状態方程式と観測方程式)の特徴についてまとめておく。
【0252】
<発明手法1>
状態方程式は、時刻n+1のK個の原画像(画素情報)が、時刻nのK個の原画像(画素情報)と、m個の原画像(画素情報)を含む有色性駆動源とから構成されている。また、観測方程式は、時刻nのK個の原画像(画素情報)が、ぼけと雑音の影響を受けてK個の劣化画像(画素情報)から構成されている。このような状態空間モデルにより、発明手法1では、観測方程式のK個の劣化画像(画素情報)を用いて、状態方程式の時刻n+1におけるK個の原画像(画素情報)を予測している。
【0253】
<発明手法2>
発明手法1において、必要な復元画像(画素情報)は1個である。それゆえ、発明手法2では、演算量を軽減するために、発明手法1の観測方程式の1個の劣化画像(画素情報)のみに着目して、発明手法1の観測方程式を再構成する。ただし、発明手法2の状態方程式は、発明手法1と同じである。このような状態空間モデルにより、発明手法2の観測方程式では、時刻nの1個の劣化画像(画素情報)のみを用いるため、過去の劣化画像(画素情報)に依存することなく、最も少ない演算量で1個の原画像(画素情報)を予測している。
【0254】
<発明手法3>
少ない演算量で高性能な画像復元能力を得るために、発明手法3の観測方程式は、有色性駆動源の影響を含む時刻nの原画像(画素情報)からなる発明手法1の観測方程式の一部と、発明手法2の観測方程式とから構成されている。このような状態空間モデルにより、発明手法3では、有色性駆動源の影響を最大限に考慮しているため、演算量の軽減を図りながら高性能な画像復元を実現している。
【0255】
<発明手法4>
処理効率を向上させるために、2次元に拡張した状態方程式は、時刻n+1のK個の原画像(画素情報)が、時刻nのK個の原画像(画素情報)と、m×q個の原画像(画素情報)を含む有色性駆動源とから構成されている。また、2次元に拡張した観測方程式は、時刻nのK個の原画像(画素情報)が、ぼけと雑音の影響を受けてK個の劣化画像(画素情報)から構成されている。このような状態空間モデルにより、発明手法4では、観測方程式のK個の劣化画像(画素情報)を用いて、状態方程式の時刻n+1におけるK個の原画像(画素情報)を予測している。
【0256】
<発明手法5>
2次元に拡張した発明手法4において、必要な復元画像(画素情報)はp×q個である。それゆえ、発明手法2では、演算量を軽減するために、発明手法4の観測方程式のp×q個の劣化画像(画素情報)のみに着目して、発明手法4の観測方程式を再構成する。ただし、発明手法5の状態方程式は、発明手法4と同じである。このような状態空間モデルにより、発明手法5の観測方程式では、時刻nのp×q個の劣化画像(画素情報)のみを用いるため、過去の劣化画像(画素情報)に依存することなく、最も少ない演算量でp×q個の原画像(画素情報)を予測している。
【0257】
<発明手法6>
少ない演算量で高性能な画像復元能力を得るために、発明手法6の観測方程式は、有色性駆動源の影響を含む時刻nの原画像(画素情報)からなる発明手法4の観測方程式の一部と、発明手法5の観測方程式とから構成されている。このような状態空間モデルにより、発明手法6では、有色性駆動源の影響を最大限に考慮しているため、演算量を軽減しながら高性能な画像復元を実現している。
【0258】
このように、発明手法1〜6によれば、いずれも、AR係数を必要としない新しい状態空間モデル(状態方程式と観測方程式)を構成して1段階処理の新しい予測法によって画像復元を行うため、1枚の劣化画像のみを用いてぼけと雑音を同時に取り除くことができる。
【0259】
また、発明手法1〜6によれば、いずれも、AR係数を必要としないため、演算量を軽減することができる。
【0260】
また、発明手法1〜3を2次元に拡張した発明手法4〜6によれば、注目画素領域および周辺画素領域のサイズをそれぞれ任意に変化させることができるため、例えば、劣化画像の性質によって1次元および2次元の新しい予測法に基づく画像復元手法を使い分けることが可能となり、効率のよい画像復元を実現することができる。
【0261】
特に、発明手法4〜6は、2次元に拡張された新しい予測法を用いた劣化画像復元手法である。この発明手法4〜6の最大の特徴は、劣化画像の性質やユーザの目的などに応じて、発明手法1、発明手法1の演算量を軽減した発明手法2、および性能を優先させた発明手法3を実現可能となるだけでなく、2次元に拡張した発明手法4、発明手法4の演算量を軽減した発明手法5、および性能を優先させた発明手法6と適応的に変更できることである。
【0262】
(実施の形態2)
実施の形態2は、状況に応じて複数の画像復元手法を適応的に切り替えて実行させる場合である。
【0263】
図115は、本発明の実施の形態2に係る画像復元装置の構成を示すブロック図である。なお、この画像復元装置200は、図1に示す画像復元装置100と同様の基本的構成を有しており、同一の構成要素にば同一の符号を付し、その説明を省略する。
【0264】
図115に示す画像復元装置200は、複数の画像復元アルゴリズムを内蔵しており、ユーザの選択によって、または、状況に応じて自動的に、内蔵した画像復元アルゴリズムを切り替えて実行可能である。複数の画像復元アルゴリズムとしては、例えば、上記した発明手法1〜6の全部または一部を用いることができる。発明手法1〜6は、上記のように、復元能力や処理速度(演算量)、処理効率、柔軟性などの点で、それぞれ特徴を有するため、適用する発明手法を状況に応じて適宜選択することによって、状況に応じた最適な画像復元を行うことができる。なお、内蔵する複数の画像復元アルゴリズムは、発明手法1〜6に限定されるわけではなく、発明手法1〜6と従来手法(例えば、上記の従来手法1〜4)との任意の組み合わせであってもよい。
【0265】
このような機能を実現するため、画像復元装置200は、復元モード指定部210と画像復元処理部220を有する。復元モード指定部210は、例えば、操作部130に設けられている。復元モード指定部210は、使用者の入力操作により、後述する復元モードを指定する。また、画像復元処理部220は、例えば、パラメータ設定部222、被写体領域判定部224を有する。また、記憶部150において、例えば、主記憶装置152は、取込画像メモリ230、被写体領域メモリ232、処理パラメータメモリ234、および復元画像メモリ236を有し、表示メモリ156は、取込画像メモリ240、処理内容メモリ242、および復元画像メモリ244を有する。
【0266】
図116は、復元モード指定部210によって指定可能な復元モードの具体例を示す図である。
【0267】
復元モード指定部210によって指定可能な復元モードは、図116に示すように、復元モード指定部210の内部構造として与えられる。図116に示す例では、復元モードとして、まず大別して、例えば、「カメラモード」、「取込画像」、「品質」が設けられている。そして、カメラモードには、例えば、「スポーツ」や「風景」、「夜景」、「接写」、「ポートレー」などのモードが設けられ、取込画像には、例えば、「写真」や「文字」、「風景」などのモードが設けられている。また、品質には、「画質優先」や「オート」、「処理速度優先」などのモードが設けられている。画質優先のモードは、さらに、「非常に細かい」、「細かい」、「普通」の3段階に分かれている。それぞれの復元モードには、予め、最適な画像復元アルゴリズムが対応付けられている。例えば、「品質」の「処理速度優先」には、発明手法5が対応付けられ、「品質」の「非常に細かい」には、発明手法3または発明手法6が対応付けられている。また、復元モードは、異なる項目から複数同時に指定することもできる。例えば、「取込画像」から「写真」を指定し、「品質」から「画質優先」の「非常に細かい」を指定した場合には、これらの条件にマッチした復元処理を行う画像復元アルゴリズムが実行される。
【0268】
画像復元処理部220において、パラメータ設定部222は、操作部130のパラメータ設定部132によって入力された各種パラメータの値を内部的に設定する。また、被写体領域判定部224は、所定の被写体領域判定アルゴリズムによって、被写体の領域を判定する。ここでは、被写体領域判定アルゴリズムは、特に限定されず、周知の任意のアルゴリズムを用いることができる。
【0269】
また、主記憶装置152において、取込画像メモリ230は、画像入力装置110によって取り込まれた画像を格納し、被写体領域メモリ232は、被写体領域判定部224によって判定された被写体の領域を格納する。被写体領域メモリ232に格納される被写体領域は、被写体の画像そのものでも被写体の画像の格納場所でもよい。処理パラメータメモリ234は、操作部130によって設定された復元処理の内容(設定されたパラメータや復元モードなど)を格納する。復元画像メモリ234は、画像復元処理部220によって復元された画像を格納する。また、表示メモリ156において、取込画像メモリ240は、主記憶装置152の取込画像メモリ230に格納された取込画像のうち、表示する取込画像を一時格納し、処理内容メモリ242は、処理パラメータメモリ234に格納された復元処理の内容(設定されたパラメータや復元モードなど)のうち、表示する処理内容を一時格納する。復元画像メモリ244は、主記憶装置152の復元画像メモリ234に格納された復元画像のうち、表示する復元画像を一時格納する。
【0270】
次いで、上記構成を有する画像復元装置200の動作について、図117に示すフローチャートを用いて説明する。なお、図117に示すフローチャートは、記憶部150(主記憶装置152または補助記憶装置154)に制御プログラムとして格納されており、図示しないCPUによって実行される。
【0271】
まず、ステップS8000では、画像入力装置110から復元処理対象の画像データ(劣化画像)を取り込み、主記憶装置152の取込画像メモリ230に格納する。
【0272】
そして、ステップS8050では、ステップS8000で取り込んだ画像データを、少なくとも1フレーム、一時的に記憶する。ここでは、取り込んだ画像データ(少なくとも1フレーム)を、1フレーム以上の画像データを一時的に蓄積するメモリ(例えば、主記憶装置152の図示しない作業用メモリ)に記憶する。
【0273】
そして、ステップS8100では、画像復元処理部220の被写体領域判定部224で、ステップS8050で一時的に記憶した取込画像データに対して、所定のレベルでプレスキャン(被写体認識)する。具体的には、取込画像データ中で所定レベル以上の変化点抽出する。
【0274】
そして、ステップS8150では、ステップS8100の結果を用いて、所定レベル以上の変化点を抽出した位置を、取込画像データに対応付けて、主記憶装置152の被写体領域メモリ232に記憶する。
【0275】
そして、ステップS8200では、プレスキャンを終了するか否かを判断する。この判断は、プレスキャンを行っていないフレームが存在するか否かに基づいて行われる。すなわち、複数のフレームが存在する場合には、フレームの数だけプレスキャンを繰り返すことになる。この判断の結果として、プレスキャンを終了する場合は(S8200:YES)、ステップS8250に進み、プレスキャンを終了しない場合は(S8200:NO)、ステップS8100に戻る。
【0276】
ステップS8250では、操作部130の復元モード指定部210および領域指定部134で、復元モードの設定と境界領域の設定を行う。設定可能な復元モードは、例えば、上記した図116に示す通りである。また、境界領域の設定では、ステップS8100でプレスキャンしたときの境界領域の確認を行う。
【0277】
そして、ステップS8300では、被写体領域か否かの判断を行う。具体的には、処理対象ブロック(注目画素とその周辺画素)が被写体領域を含んでいるか否かを判断する。この判断は、例えば、プレスキャンした段階で、劣化画像の中のどの位置にユーザが注目している被写体があるのかを画素の輝度値により自動的に判定しておき、この結果を利用したり、あるいは、ユーザにタッチパネルを使って選択的に被写体領域を判定させたりすることによって行われる。この判断の結果として、処理対象ブロックが被写体領域を含む場合は(S8300:YES)、ステップS8350に進み、処理対象ブロックが被写体領域を含まない場合は(S8300:NO)、ステップS8400に進む。
【0278】
ステップS8350では、実行する画像復元手法として、ステップS8250で設定した復元モードに対応する1次元の画像復元手法(例えば、発明手法1〜3のいずれか)を選択する。また、ステップS8400では、実行する画像復元手法として、ステップS8250で設定した復元モードに対応する2次元の画像復元手法(例えば、発明手法4〜6のいずれか)を選択する。すなわち、本実施の形態では、処理対象ブロックが被写体領域を含むか否かで大きく1次元復元手法を実行するか2次元復元手法を実行するかを選択し、その後、選択された1次元復元手法または2次元復元手法の中から、例えば、取込画像(写真や文字など)、品質(画質やスピードなど)、または様々な用途(ポートレートや夜景など)に適応した復元モードに応じて、最適な復元手法を選択するようにしている。
【0279】
そして、ステップS8450では、画像復元処理部220で、ステップS8350またはステップS8400で選択した画像復元手法を用いて、設定されたパラメータに従った画像復元処理を行う。ここで、パラメータの設定は、操作部130のパラメータ設定部132で行われる。
【0280】
そして、ステップS8500では、当該フレームについて画像復元処理が終了したか否かを判断する。この判断の結果として、当該フレームについて画像復元処理が終了した場合は(S8500:YES)、ステップS8550に進み、当該フレームについて画像復元処理が終了していない場合は(S8500:NO)、ステップS8300に戻る。
【0281】
ステップS8550では、全フレームについて画像復元処理が終了したか否かを判断する。この判断の結果として、全フレームについて画像復元処理が終了していない場合は(S8550:NO)、ステップS8250に戻り、全フレームについて画像復元処理が終了した場合は(S8550:YES)、上記一連の処理を終了する。したがって、複数のフレームが存在する場合には、フレームの数だけ画像復元処理を繰り返すことになる。
【0282】
このように、本実施の形態によれば、状況に応じて複数の画像復元手法を適応的に切り替えて実行させるため、復元する画素に適した画像復元処理を、1チップで、高精度かつ高速に行うことができる。
【産業上の利用可能性】
【0283】
本発明に係る画像復元装置および画像復元方法は、画像復元能力を向上することができるシンプルで実用的な画像復元装置および画像復元方法として有用である。
【符号の説明】
【0284】
100、200 画像復元装置
110 画像入力装置
112 カメラ
114 スキャナ
116、186 記録メディア
118、188 モデム
120 入力インタフェース部
130 操作部
132、222 パラメータ設定部
134 領域指定部
136、210 復元モード指定部
140 内部インタフェース部
150 記憶部
152 主記憶装置
154 補助記憶装置
156 表示メモリ
160、220 画像復元処理部
160a 第1復元処理部
160b 第2復元処理部
162 初期設定部
164 相関演算部
166 重み係数算出部
168 最適推定値算出部
170 出力インタフェース部
180 画像出力装置
182 ディスプレイ
184 プリンタ
224 被写体領域判定部
230、240 取込画像メモリ
232 被写体領域メモリ
234 処理パラメータメモリ
236、244 復元画像メモリ
242 処理内容メモリ

【特許請求の範囲】
【請求項1】
原画像情報に不必要な情報が混在した劣化画像情報のみから前記原画像情報を推定する画像復元装置であって、
時刻nのみの劣化画像情報に対して、時刻nまたは時刻n+1までの情報により前記原画像情報を含む時刻n+1でのシステムの状態量を推定した場合の推定誤差の相関値を算出する相関演算部と、
時刻nのみの劣化画像情報に対して、前記相関演算部によって算出された相関値を用いて、時刻n+1までの情報による時刻n+1での前記状態量の最適推定値と、時刻nまでの情報による時刻n+1での前記状態量の最適推定値と、前記劣化画像情報を含む観測量の推定誤差と、の関係を規定するための重み係数を算出する重み係数算出部と、
時刻nのみの劣化画像情報に対して、前記重み係数算出部によって算出された重み係数を用いて、時刻nまたは時刻n+1までの情報による時刻n+1での前記状態量の最適推定値を算出する最適推定値算出部と、
を有する画像復元装置。
【請求項2】
原画像情報に不必要な情報が混在した劣化画像情報のみから前記原画像情報を推定する画像復元方法であって、
時刻nのみの劣化画像情報に対して、時刻nまたは時刻n+1までの情報により前記原画像情報を含む時刻n+1でのシステムの状態量を推定した場合の推定誤差の相関値を算出する相関演算工程と、
時刻nのみの劣化画像情報に対して、前記相関演算工程で算出した相関値を用いて、時刻n+1までの情報による時刻n+1での前記状態量の最適推定値と、時刻nまでの情報による時刻n+1での前記状態量の最適推定値と、前記劣化画像情報を含む観測量の推定誤差と、の関係を規定するための重み係数を算出する重み係数算出工程と、
時刻nのみの劣化画像情報に対して、前記重み係数算出工程で算出した重み係数を用いて、時刻nまたは時刻n+1までの情報による時刻n+1での前記状態量の最適推定値を算出する最適推定値算出工程と、
を有する画像復元方法。
【請求項3】
原画像情報に不必要な情報が混在した劣化画像情報のみから前記原画像情報を推定するための画像復元プログラムであって、
コンピュータに、
時刻nのみの劣化画像情報に対して、時刻nまたは時刻n+1までの情報により前記原画像情報を含む時刻n+1でのシステムの状態量を推定した場合の推定誤差の相関値を算出する相関演算工程と、
時刻nのみの劣化画像情報に対して、前記相関演算工程で算出した相関値を用いて、時刻n+1までの情報による時刻n+1での前記状態量の最適推定値と、時刻nまでの情報による時刻n+1での前記状態量の最適推定値と、前記劣化画像情報を含む観測量の推定誤差と、の関係を規定するための重み係数を算出する重み係数算出工程と、
時刻nのみの劣化画像情報に対して、前記重み係数算出工程で算出した重み係数を用いて、時刻nまたは時刻n+1までの情報による時刻n+1での前記状態量の最適推定値を算出する最適推定値算出工程と、
を実行させるための画像復元プログラム。

【図1】
image rotate

【図2】
image rotate

【図3】
image rotate

【図4】
image rotate

【図5】
image rotate

【図6】
image rotate

【図7】
image rotate

【図8】
image rotate

【図9】
image rotate

【図10】
image rotate

【図11】
image rotate

【図12】
image rotate

【図13】
image rotate

【図14】
image rotate

【図15】
image rotate

【図16】
image rotate

【図17】
image rotate

【図18】
image rotate

【図19】
image rotate

【図20】
image rotate

【図21】
image rotate

【図22】
image rotate

【図23】
image rotate

【図24】
image rotate

【図25】
image rotate

【図26】
image rotate

【図27】
image rotate

【図28】
image rotate

【図29】
image rotate

【図30】
image rotate

【図31】
image rotate

【図32】
image rotate

【図33】
image rotate

【図34】
image rotate

【図35】
image rotate

【図36】
image rotate

【図37】
image rotate

【図38】
image rotate

【図39】
image rotate

【図40】
image rotate

【図41】
image rotate

【図42】
image rotate

【図43】
image rotate

【図44】
image rotate

【図45】
image rotate

【図46】
image rotate

【図47】
image rotate

【図48】
image rotate

【図49】
image rotate

【図50】
image rotate

【図51】
image rotate

【図52】
image rotate

【図53】
image rotate

【図54】
image rotate

【図55】
image rotate

【図56】
image rotate

【図57】
image rotate

【図58】
image rotate

【図59】
image rotate

【図60】
image rotate

【図61】
image rotate

【図62】
image rotate

【図63】
image rotate

【図64】
image rotate

【図65】
image rotate

【図66】
image rotate

【図67】
image rotate

【図68】
image rotate

【図69】
image rotate

【図70】
image rotate

【図71】
image rotate

【図72】
image rotate

【図73】
image rotate

【図74】
image rotate

【図75】
image rotate

【図76】
image rotate

【図77】
image rotate

【図78】
image rotate

【図79】
image rotate

【図80】
image rotate

【図81】
image rotate

【図82】
image rotate

【図83】
image rotate

【図84】
image rotate

【図85】
image rotate

【図86】
image rotate

【図87】
image rotate

【図88】
image rotate

【図89】
image rotate

【図90】
image rotate

【図91】
image rotate

【図92】
image rotate

【図93】
image rotate

【図94】
image rotate

【図95】
image rotate

【図96】
image rotate

【図97】
image rotate

【図98】
image rotate

【図99】
image rotate

【図100】
image rotate

【図101】
image rotate

【図102】
image rotate

【図103】
image rotate

【図104】
image rotate

【図105】
image rotate

【図106】
image rotate

【図107】
image rotate

【図108】
image rotate

【図109】
image rotate

【図110】
image rotate

【図111】
image rotate

【図112】
image rotate

【図113】
image rotate

【図114】
image rotate

【図115】
image rotate

【図116】
image rotate

【図117】
image rotate


【公開番号】特開2011−165013(P2011−165013A)
【公開日】平成23年8月25日(2011.8.25)
【国際特許分類】
【出願番号】特願2010−28278(P2010−28278)
【出願日】平成22年2月10日(2010.2.10)
【出願人】(803000115)学校法人東京理科大学 (545)
【Fターム(参考)】