(57)【要約】  同一の対象を撮影した複数枚の観測画像をもとに、より高い解像度の画像を生成する超解像問題において、隠れマルコフ確率場事前分布とその周辺化をとりこんだアルゴリズムにより高解像度画像を推定する方法およびプログラムを提供する。劣化された複数の観測画像から原画像を復元する解像において、複数の観測画像についてベイズ推定の手法を用い、特に、画像の生成過程に、局所的な滑らかさと大域的な滑らかさの２つからなる階層的な事前分布を仮定する手法により、ボケた画像を推定してしまうという問題を回避し、またパラメータの推定精度を向上する。また、エッジ保存性を持つ隠れマルコフ確率場を事前分布とし、かつ周辺化をおこなうことで高精度の劣化過程の同定と高解像度の画像推定の両方を実現する。さらに、適切な近似により多項式オーダとして計算コストを小さくする。
2007122838.tif 000002
【特許請求の範囲】
【請求項１】
原画像に関連する複数の低解像度の画像から高解像度の画像を生成するための方法であって、
（１）前記複数の低解像度の画像を入力するステップと、
（２）前記高解像度の画像の事前分布を局所的な滑らかさ情報（各画素は近傍の画素と相関を有する）および大域的な平滑化情報（各画素は近傍の画素と全く相関が無い）とを組み合わせた階層化分布と定義するステップと、
（３）前記高解像度の画像の事後分布をベイズ推定の手法を用いて推定するステップと、
（４）原画像の劣化過程のパラメータをベイズ推定の手法を用いて推定するステップと、
（５）ステップ（３）で推定した高解像度の画像の事後分布に基づいて、ステップ（４）で収束したパラメータを用い、高解像度の画像を生成するステップとを含むことを特徴とする超解像法。
【請求項２】
原画像に関連する複数の低解像度の画像から高解像度の画像を生成するための方法であって、
（１）前記複数の低解像度の画像を入力するステップと、
（２）前記高解像度の画像の事前分布を局所的な滑らかさ情報（各画素は近傍の画素と相関を有する）および大域的な平滑化情報（各画素は近傍の画素と全く相関が無い）とを組み合わせた階層化分布と定義するステップと、
（３）前記高解像度の画像の事後分布をベイズ推定の手法を用いて推定するステップと、
（４）原画像の劣化過程のパラメータをベイズ推定の手法を用いて推定するステップと、
（５）前記パラメータの値が収束するまで、ステップ（３）〜（４）を反復するステップと、
（６）ステップ（３）で推定した高解像度の画像の事後分布に基づいて、ステップ（４）で収束したパラメータを用い、高解像度の画像を生成するステップとを含むことを特徴とする超解像法。
【請求項３】
ステップ（２）における高解像度の画像の事前分布を、エッジ保存性を有する隠れマルコフ確率場と定義することを特徴とする請求項１又は２に記載の超解像法。
【請求項４】
前記ベイズ推定の手法を用いて推定するステップ（３）およびステップ（４）において、推定演算が多項式オーダとなるように所定の近似を行うことを特徴とする請求項３に記載の超解像法。
【請求項５】
コンピュータに、原画像に関連する複数の低解像度の画像から高解像度の画像を生成させるためのプログラムであって、
（１）原画像の複数の低解像度の画像を入力させる手順と、
（２）高解像度の画像の事前分布を局所的な滑らかさ情報（各画素は近傍の画素と相関を有する）および大域的な平滑化情報（各画素は近傍の画素と全く相関が無い）とを組み合わせた階層化分布と定義させる手順と、
（３）高解像度の画像の事後分布をベイズ推定アルゴリズムにより推定させる手順と、
（４）原画像の劣化過程のパラメータをベイズ推定アルゴリズムにより推定させる手順と、
（５）手順（３）で推定した高解像度の画像の事後分布に基づいて、手順（４）で収束したパラメータを用い、高解像度の画像を生成させる手順を、
を実行させるための超解像プログラム。
【請求項６】
コンピュータに、原画像に関連する複数の低解像度の画像から高解像度の画像を生成させるためのプログラムであって、
（１）原画像の複数の低解像度の画像を入力させる手順と、
（２）高解像度の画像の事前分布を局所的な滑らかさ情報（各画素は近傍の画素と相関を有する）および大域的な平滑化情報（各画素は近傍の画素と全く相関が無い）とを組み合わせた階層化分布として判断させる手順と、
（３）高解像度の画像の事後分布をベイズ推定アルゴリズムにより推定させる手順と、
（４）原画像の劣化過程のパラメータをベイズ推定アルゴリズムにより推定させる手順と、
（５）前記パラメータの値が収束するまで、手順（３）〜（４）を反復させる手順と、
（６）手順（３）で推定した高解像度の画像の事後分布に基づいて、手順（４）で収束したパラメータを用い、高解像度の画像を生成させる手順と、
を実行させるための超解像プログラム。
【請求項７】
手順（２）における高解像度の画像の事前分布を、エッジ保存性を有する隠れマルコフ確率場と定義させることを特徴とする請求項５又は６に記載の超解像プログラム。
【請求項８】
前記ベイズ推定アルゴリズムにより推定する手順（３）および手順（４）において、推定演算が多項式オーダとなるように所定の近似演算を行わせる手順を実行させるための請求項７に記載の超解像プログラム。
【発明の詳細な説明】
【技術分野】
【０００１】
本発明は、画像処理に関し、特に複数の低解像度の画像を使用して劣化画像から原画像を高精度に復元する超解像の方法およびプログラムに関するものである。
【背景技術】
【０００２】
超解像とは，カメラ、ビデオなどレンズやＣＣＤ（Charge Coupled Devices）で構成される結像系に内在する物理的限界を超えて，画像の解像度を高めることをいう。カメラ、ビデオなど撮像装置では、レンズやＣＣＤなどの構成部品の品質が撮影される画像の解像度を決定する。また、レンズの収差に起因するボケやＣＣＤの密度に依存するサブサンプリング，さらにＣＣＤのショットノイズなど観測に混入するノイズにより，撮影された画像は、真の撮影対象の画像（原画像）よりも劣化したものになる．【０００３】
レンズやＣＣＤなどの構成部品のハードウエアの工夫では超えようのない物理的な限界を突破した高解像度の画像を得るために、ソフトウエア的な超解像処理を行い、カメラ、ビデオなどの撮像装置の能力を超える高解像度の画像を得る方法が研究されている。また、撮像対象の物体を僅かずつ異なる位置・角度から撮影した画像を複数枚用意できる場合、それら複数の観測画像の情報を統計的に情報処理することで、より高い解像度の画像を生成できることが知られている。【０００４】
低解像度の画像を使用して劣化画像から原画像を高精度に復元する超解像手法のうち統計的な解法に関する従来の研究は、最尤法，ＭＡＰ（maximum a posteriori）法，ベイズ統計法の３つに大別される。
超解像のアリゴリズムは数多くあるが、ほとんどは位置ズレ推定と高解像度画像推定は別のモデルにもとづいている。同一の確率モデルによる超解像法では、ＭＡＰ法を用いたHardie らや（非特許文献１）、ベイズ統計法を採用したTipping＆Bishopの方法が知られている（非特許文献２）。非特許文献２では、高解像度画像の周辺化を用いたベイズ統計法により、モデルパラメタの推定精度が著しく上昇することが示されている。Tipping＆Bishopは、高解像度画像の事前分布に単純なガウス確率場を採用し、原画像の局所なめらかさ制約を事前分布として加え、事後確率最大化推定を行い、一般には困難な周辺化演算を解析的に実行可能にしている。【０００５】
【非特許文献１】R. C. Hardie, K. J. Barnard, and E. E. Armstrong, “Joint MAP registration and highresolution image estimation using a sequence of undersampled images,” IEEE Transactions on Image Processing, vol. 6, no. 12, pp. 16211633, 1997.
【非特許文献２】M. E. Tipping and C. M. Bishop, “Bayesian image super-resolution,” in Advances in Neural Processing Systems 15 (S. Becker, S. Thrun, and K. Obermayer, eds.), pp. 12791286, MIT Press, 2003.
【発明の開示】
【発明が解決しようとする課題】
【０００６】
Tipping＆Bishopは、高解像度画像の周辺化を用いたベイズ統計法を超解像度に最初に適用し、周辺化によりモデルパラメタの推定精度が飛躍的に向上することを示した。Tipping＆Bishopは、ガウス確率場と呼ばれる確率過程を原画像に仮定した（高解像度画像の事前分布に用いた単純なガウス確率場を用いた）ものであり、これは局所的な滑らかさ制約を加えたことに相当する。しかし、この滑らかさに対する強い制約のため、高解像度画像の推定の際にはボケた画像を推定してしまうという問題があった。
また、原画像の真の事前分布から離れた単純なガウス確率場の事前分布を用いてモデルパラメタの推定を行うため、モデルパラメタの推定精度が悪くなるという問題があった。【０００７】
現実の画像は物体の境界にエッジが存在しており、局所的な滑らかさばかりでなく、近傍の画素とは相関のない大域的な滑らかさもあるので、両者を制御する隠れ変数を持ち、両者の特性を併せ持つ階層的な事前分布がより適切と考えられるが、こうした階層的な事前分布を用いると推定（画像復元）に多大な計算量が必要となる問題点があった。【０００８】
上記問題点に鑑み、本発明は、現実の画像の性質に関する事前分布として、局所的な滑らかさと大域的な滑らかさの２つからなる階層的な事前分布を用いて、また、原画像の劣化過程におけるパラメータ推定にベイズ推定を用いることで、高精度の超解像を実現することを目的とする。
【課題を解決するための手段】
【０００９】
上記目的を達成するため、本発明の超解像法は、原画像に関連する複数の低解像度の画像から高解像度の画像を生成するための方法であって、
（１）前記複数の低解像度の画像を入力するステップと、
（２）前記高解像度の画像の事前分布を局所的な滑らかさ情報（各画素は近傍の画素と相関を有する）および大域的な平滑化情報（各画素は近傍の画素と相関をもたない）とを組み合わせた階層化分布と定義するステップと、
（３）前記高解像度の画像の事後分布をベイズ推定の手法を用いて推定するステップと、
（４）原画像の劣化過程のパラメータをベイズ推定の手法を用いて推定するステップと、
（５）ステップ（３）で推定した高解像度の画像の事後分布に基づいて、ステップ（４）で収束したパラメータを用い、高解像度の画像を生成するステップと、を含むことを特徴とするものである。【００１０】
原画像に対して、位置ずれ、回転、ボケ、低解像化、ノイズなどが加わることで、観測画像が得られるというのは一般的な状況である。このように劣化された観測画像から原画像を復元する解像は、劣化が不可逆過程であるため、たとえ劣化過程が既知であったとしても不可能であるが、観測画像が複数あり劣化過程が観測画像ごとに独立に近い場合は、それら複数の観測画像について一種の多数決を用いることで推定可能となる。【００１１】
本発明の超解像法は、この多数決に統計的推定法の一種であるベイズ推定の手法を用いることで、高精度の解像が可能となる超解像を実現する。特に、画像の生成過程に、局所的な滑らかさと大域的な滑らかさの２つからなる階層的な事前分布を仮定する手法により、高解像度画像の推定の際にはボケた画像を推定してしまうという問題を回避する。
また、真の事前分布に近い、局所的な滑らかさと大域的な滑らかさの２つからなる階層的な事前分布を仮定することで、パラメータの推定精度を向上する。
本発明の超解像法では、原画像の劣化過程のパラメータ（位置ズレなどのパラメータ）については陽に推定し、位置ズレなどのパラメータを推定する際に、高解像度画像に関しては陽に定めず確率分布を推定し、それらの期待値をとるような周辺化によってオーバーフィッティングを防いでいるが、位置ズレなどのパラメータに関しても陽に定めず、確率分布を推定するようなベイズ統計法を用いることも可能である。【００１２】
ここで、局所的な滑らかさと大域的な滑らかさの２つからなる階層的な事前分布とは、原画像の境界を表現し得る混合分布のことで、各画素は近傍の画素と相関を有するか、或いは、全く違う値をとってもよいとする原画像の性質に関する知識である。【００１３】
また、本発明の超解像法は、原画像に関連する複数の低解像度の画像から高解像度の画像を生成するための方法であって、
（１）前記複数の低解像度の画像を入力するステップと、
（２）前記高解像度の画像の事前分布を局所的な滑らかさ情報（各画素は近傍の画素と相関を有する）および大域的な平滑化情報（各画素は近傍の画素と相関をもたない）とを組み合わせた階層化分布と定義するステップと、
（３）前記高解像度の画像の事後分布をベイズ推定の手法を用いて推定するステップと、
（４）原画像の劣化過程のパラメータをベイズ推定の手法を用いて推定するステップと、
（５）前記パラメータの値が収束するまで、ステップ（３）〜（４）を反復するステップと、
（６）ステップ（３）で推定した高解像度の画像の事後分布に基づいて、ステップ（４）で収束したパラメータを用い、高解像度の画像を生成するステップと、を含むことを特徴とするものである。【００１４】
本発明においては、高解像度画像を隠れ変数とし、その確率分布(のパラメータ)と、位置ズレなどのパラメータと交互に最適化する方法(ＥＭアルゴリズム)を用いているが、これら2種類のパラメータは同時に推定することも、交互に最適化することも可能である。どちらの場合でもこれら2種類のパラメータに関しての自由エネルギの最小化となるようにパラメータを収束させるのである。【００１５】
また、本発明の超解像法は、上記のステップ（２）における高解像度の画像の事前分布を、エッジ保存性を有する隠れマルコフ確率場と定義することを特徴とする。
エッジ保存性を持つ隠れマルコフ確率場を事前分布とし、かつ周辺化をおこなうことで高精度の劣化過程の同定と高解像度の画像推定の両方を実現するのである。【００１６】
また、本発明の超解像法は、上記のベイズ推定の手法を用いて推定するステップ（３）およびステップ（４）において、推定演算が多項式オーダとなるように所定の近似を行うことを特徴とする。
階層的な事前分布を仮定することで、推定が指数オーダとなり推定に多大な計算量が必要となるが、適切な近似により多項式オーダとして計算コストを小さくするのである。
ここで、適切な近似とは、事後分布を画素ごとに完全に独立なmixture of Gaussiansで置き換えるという近似手法である（詳細については後述する）。なお、画素ごとに完全に独立なmixture of Gaussiansで置き換えなくとも、近傍画素間だけに相関を許すといった近似計算も可能である。【００１７】
次に、本発明の超解像プログラムは、コンピュータに、原画像に関連する複数の低解像度の画像から高解像度の画像を生成させるためのプログラムであって、
（１）原画像の複数の低解像度の画像を入力させる手順と、
（２）高解像度の画像の事前分布を局所的な滑らかさ情報（各画素は近傍の画素と相関を有する）および大域的な平滑化情報（各画素は近傍の画素と全く相関が無い）とを組み合わせた階層化分布と定義させる手順と、
（３）高解像度の画像の事後分布をベイズ推定アルゴリズムにより推定させる手順と、
（４）原画像の劣化過程のパラメータをベイズ推定アルゴリズムにより推定させる手順と、
（５）手順（３）で推定した高解像度の画像の事後分布に基づいて、手順（４）で収束したパラメータを用い、高解像度の画像を生成させる手順を、を実行させることを特徴とする。【００１８】
また、本発明の超解像プログラムは、コンピュータに、原画像に関連する複数の低解像度の画像から高解像度の画像を生成させるためのプログラムであって、
（１）原画像の複数の低解像度の画像を入力させる手順と、
（２）高解像度の画像の事前分布を局所的な滑らかさ情報（各画素は近傍の画素と相関を有する）および大域的な平滑化情報（各画素は近傍の画素と全く相関が無い）とを組み合わせた階層化分布として判断させる手順と、
（３）高解像度の画像の事後分布をベイズ推定アルゴリズムにより推定させる手順と、
（４）原画像の劣化過程のパラメータをベイズ推定アルゴリズムにより推定させる手順と、
（５）前記パラメータの値が収束するまで、手順（３）〜（４）を反復させる手順と、
（６）手順（３）で推定した高解像度の画像の事後分布に基づいて、手順（４）で収束したパラメータを用い、高解像度の画像を生成させる手順と、を実行させることを特徴とする。【００１９】
また、上記の超解像プログラムの手順（２）における高解像度の画像の事前分布を、エッジ保存性を有する隠れマルコフ確率場と定義することが好ましい。エッジ保存性を持つ隠れマルコフ確率場を事前分布とし、かつ周辺化をおこなうことで高精度の劣化過程の同定と高解像度の画像推定の両方を実現するためである。【００２０】
さらに、上記の超解像プログラムのベイズ推定アルゴリズムにより推定する手順（３）および手順（４）において、推定演算が多項式オーダとなるように所定の近似演算を行う手順を実行させることが好ましい。ベイズ推定アルゴリズムによる推定が指数オーダとなり推定に多大な計算量が必要となるところを、適切な近似により多項式オーダとして計算コストを小さくすることができるためである。
【発明の効果】
【００２１】
本発明の超解像法および超解像プログラムによれば、輪郭がはっきりとした高解像度画像を得ることができる。
また、本発明の超解像法および超解像プログラムによれば、全ての計算を多項式オーダに収める近似手法を用いることで、計算コストを減少することができる。
【発明を実施するための最良の形態】
【００２２】
以下、図面を参照しつつ、本発明の超解像法および超解像プログラムの実施例を説明する。図１に、本発明の超解像法のフローチャートを示す。【００２３】
（１）先ず、２個以上の低解像度の画像を入力する（ステップＳ１０）。ここで、低解像度の画像とは、原画像に対して、位置ずれ、回転、ボケ、ノイズなどが印加されたり、レンズの収差に起因するボケやＣＣＤの密度に依存するサブサンプリング，さらにＣＣＤのショットノイズなど観測に混入するノイズによって、劣化された画像を指している。【００２４】
（２）次に、高解像度の画像の事前分布について、局所的な滑らかさ情報（各画素は近傍の画素と相関を有する）および大域的な平滑化情報（各画素は近傍の画素と全く相関が無い）とを組み合わせた階層化分布と定義する（ステップＳ２０）。
具体的には、画像の事前分布を、下式で示されるエッジ保存性を有する隠れマルコフ確率場と定義する。【００２５】
【数１】
2007122838.tif 000003
【００２６】
ここで、Σηはすべてのηi に関する和Ση1・・・ΣηPo を表わすものである。また、隠れ変数ηiは、画素i のコネクティビティを表現するもので、ηi ＝１のときは近傍の画素値と近い値を、ηi＝０のときは近傍に関係ない値をとりやすくするものである。
ηの導入により、画素値の断裂を表現することが可能になり、エッジを保持した推定が期待できることになる。また、本事前分布はエッジを意識した分布であるため、遮蔽などによってエッジの生じる実際の画像の分布に適合しており、位置ズレの推定においても高精度を示すことが期待できるのである。【００２７】
ここで、ηの導入の利点について、経済性と平滑化の２つの側面を説明する。
経済性の側面では、画素間に隠れ変数を付与するのではなく，画素ごとに付与するηの導入は必要とされる近似の削減につながり、周辺化の計算を多項式時間で実行するために利便性が高いといった利点がある。
また、平滑化の側面では、1 つのηi は画素i の近傍の画素全てに影響を及ぶため、ηiを導入する方がより滑らかな画像になりやすく、真の分布の性質に適合するといった利点がある。【００２８】
あるηの実現が与えられたときの高解像度画像ｘの分布は、下式で示されるギブス分布をとる。【００２９】
【数２】
2007122838.tif 000004
【００３０】
ここで、Ｚpri は正規化項で、ρ（＞０）は逆温度で、Ｅは下式のように定義されるエネルギ関数である。【００３１】
【数３】
2007122838.tif 000005
【００３２】
但し、Ｓ（ｘ、η）とＧ（ｘ、η）は、下式で示されるものである。【００３３】
【数４】
2007122838.tif 000006
【００３４】
【数５】
2007122838.tif 000007
【００３５】
ｃ（＞０）により、Ｓ（ｘ、η）とＧ（ｘ、η）のバランスが調節されることになる。この第１項Ｓ（ｘ、η）は、ある画素の値が近傍の画素値と近いほど確率が高いことを要請している局所的な滑らかさ制約である。一方、第２項Ｇ（ｘ、η）は、近傍の画素値に関係無く、νに近いほど確率が高いという条件を課す大域的な平滑化項である。この２つのどちらをとるかをηで制御する。すなわち、ηi＝１のとき、画素i は第１項Ｓを小さくしようと近傍の画素と近い値をとろうとし、ηi ＝０のとき画素i は第２項を小さくしようと近傍の画素値によらずνに近い値をとろうとする。ηが与えられたとき、エネルギ関数Ｅは、常に非負、かつ、ｘの２次関数であるから、正定値対称行列Ａを用いて下式で表される。【００３６】
【数６】
2007122838.tif 000008
【００３７】
上式のエネルギ関数Ｅにおいて、各項の行列、変数を下式に示す。【００３８】
【数７】
2007122838.tif 000009
【００３９】
上記から、あるηの実現が与えられたときの高解像度画像の事前分布Ｐ（ｘ｜η）は、結局のところ、ガウス分布であり下式で表される。【００４０】
【数８】
2007122838.tif 000010
【００４１】
ここで、隠れ変数ηの事前分布Ｐ（η）は、各変数毎に独立に、確率ａでηi＝１，確率１−aでηi＝０をとる分布、すなわち下式で表されるベルヌイ分布に設定する。【００４２】
【数９】
2007122838.tif 000011
【００４３】
以上より、エッジ保存性を有する隠れマルコフ確率場と定義した画像の事前分布は、混合ガウス分布となり、混合数は２Ｐｏ，混合比は上式のBernoulli（η｜ａ）、成分は下式で示される。【００４４】
【数１０】
2007122838.tif 000012
【００４５】
（３）次に、高解像度の画像の事後分布および原画像の劣化過程のパラメータを推定する。この推定は、２Ｐｏ個の項からなる和を含むため、現実的な時間で計算するのは不可能である。
本発明は、現実的な時間で算出すべく、試験分布を計算の容易な形にとったＥＭ（expectation−maximization）アルゴリズムを用いている。
ＥＭアルゴリズムは、隠れ変数をもつ確率モデルにおいて最尤解を求めるための計算手順であり、自由エネルギのＥステップとＭステップによる交代最適化を反復するアルゴリズムである。隠れ変数は、直接観測できない高解像度画像ｘと離散変数ηであると仮定する。ＥＭアルゴリズムを用いると、直接に尤度を計算するのが困難な場合にでも、より計算の容易な自由エネルギを計算することで効率的に最尤解を求めることができるのである。【００４６】
ＥＭアルゴリズムは、自由エネルギＦ[ｑ，θ]の交代最小化として定式化される。自由エネルギＦ[ｑ，θ]は，確率密度関数ｑの汎関数かつモデルパラメタθの関数であり、下式のように定義される。【００４７】
【数１１】
2007122838.tif 000013
【００４８】
ＥＭアルゴリズムの実行に際しては、まず適当なモデルパラメタの初期値θ０を用意する。ＥＭアルゴリズムの第l（＞１）反復における計算は、下式で表されるＥステップとＭステップの２ステップからなる。収束に至るまで反復してＥステップとＭステップの計算をくりかえし、収束したときのθlをモデルパラメタの推定値とする。【００４９】
【数１２】
2007122838.tif 000014
【００５０】
図１におけるフローでは、高解像度の画像の事後分布をベイズ推定の手法を用いて推定するステップＳ３０と、原画像の劣化過程のパラメータをベイズ推定の手法を用いて推定するステップＳ４０と、パラメータの値が収束するまで、ステップＳ３０〜ステップＳ４０を反復するステップＳ５０となる。【００５１】
ステップＳ３０における高解像度の画像の事後分布の推定は、指数オーダの項からなる和を含むことから、現実的な時間で計算するのは困難である。そこで、本発明では、試験分布q を任意の確率分布にとるのではなく、自由エネルギの計算が可能であるような簡単な形に限定することでＥＭアルゴリズムを実行可能とするのである。すなわち、試験分布ｑ（η，ｘ）を次元ごとに独立に因子分解された下式で表される混合ガウス分布にとるのである。【００５２】
【数１３】
2007122838.tif 000015
【００５３】
上式により定義される試験分布は、各次元が独立であり、画素i に対し，ηi が１か０かに対応してガウス分布が２つ存在し、それぞれの成分の混合比はκiηi，平均はψiηi，分散はξiηiである。【００５４】
また、自由エネルギは、因子化試験分布ｑ（η，ｘ｜τ）のもとではもはや汎関数ではなく、試験分布のパラメータτとモデルパラメタθの関数Ｆ（τ，θ）となる。これにともないＥステップとＭステップは下式と表すことが可能となる。【００５５】
【数１４】
2007122838.tif 000016
【００５６】
この自由エネルギＦ（τ，θ）は，結合分布を連鎖律にしたがい展開して、下式の如く分解されることになる。【００５７】
【数１５】
2007122838.tif 000017
【００５８】
ここで、各項のＦｃ（τ）、Ｆｓ（τ）、Ｆｏ（τ，θ）、ＦＨ（τ）は、下式から算出される。【００５９】
【数１６】
2007122838.tif 000018
【００６０】
（６）図１に示されるフローで、最後に、ステップ（３）で推定した高解像度の画像の事後分布に基づいて、ステップ（４）で収束したパラメータを用い、高解像度の画像を生成する（ステップＳ６０）。
これにより、輪郭が明瞭な高精度の高解像度の画像を復元することができる。【００６１】
以下、実施例１で、人工的に生成した同じ観測画像データに対して本発明の超解像法を適用して、位置ズレの推定や高解像度画像の推定を行い、その性能をTipping＆Bishopのガウス確率場を事前分布とした方法と比較する。
【実施例１】
【００６２】
先ず、高解像度画像を用意し、これに平行移動，回転，ボケ，ダウンサンプリングを加えて１６枚の低解像度画像を人工的に作り、さらにガウスノイズを加えて観測画像データを生成した。平行移動量と回転角度は、それぞれ区間[−２，２]，[−４π／180，４π／180]（ラジアン）の一様分布にしたがってサンプリングした。ＰＳＦ（Point spread function）には、ガウスＰＳＦを適用し、パラメータγ は２に設定した。また、縮小率は縦横１／４倍とした。ノイズ条件のみを変えて４つのデータを作成した。ノイズ条件の基準は、下式に示されるＢＳＮＲ (Blurred Signal to Noise Ratio)としている。【００６３】
【数１７】
2007122838.tif 000019
【００６４】
ここで、μはノイズを加える前の低解像度劣化画像の画素値の平均値であり、下式で表される。【００６５】
【数１８】
2007122838.tif 000020
【００６６】
また、σ２はガウスノイズの分散である。３つのデータセットのノイズ条件は、ノイズ無し，ＢＳＮＲが３０ｄｂ，２０ｄｂである。【００６７】
アルゴリズムの実行に際してのパラメータの初期化は次のようにした。モデルパラメタθ＝[{St, φt}t, γ] については、Stとφtはすべてゼロに設定し、γ は４と設定した。
試験分布のパラメータτ＝[κ，ψ，ξ]については、κはすべて０．５とし、ψは１６枚の観測画像を高解像度に線型補完した画像の平均値とし、ξは区間[0.05−0.003，0.05−0.003] の一様分布から生成した。【００６８】
ここで、勾配ベースの数値最適化法には、scaled conjugate gradient 法を用いた。
実験に用いた原画像を図２に示す。これより人工的に生成された１６枚の観測画像を図３に示す。観測画像のうちの１枚を拡大した画像を図４に示す。
人工的に生成した同じ観測画像（図３に示したもの）に対し、本発明の超解像法を用いて、位置ズレの推定や高解像度画像の推定を行ったノイズ無しのデータからの超解像度画像を図５（ａ）に示す。比較のために、Tipping＆Bishopのガウス確率場を事前分布とした方法を用いて推定を行った超解像度画像を図５（ｂ）に示す。
また、本発明の超解像法を用いて推定を行った、ノイズ条件が３０ｄｂ ＢＳＮＲのデータからの超解像度画像を図６（ａ）に示す。比較のために、Tipping＆Bishopのガウス確率場を事前分布とした方法を用いて推定を行った超解像度画像を図６（ｂ）に示す。【００６９】
本発明の超解像法を用いた推定結果の超解像度画像も、Tipping＆Bishopのガウス確率場を事前分布とした方法を用いた推定結果の超解像度画像も、双方とも観察画像では読めなかった文字が充分判別できる程度までの解像度の復元が見られた。図５および図６から、本発明の超解像法を用いた推定結果の超解像度画像の方が、Tipping＆Bishopのガウス確率場を事前分布とした方法をよりシャープな画像を推定していることが理解できる。【００７０】
モデルパラメタと高解像度画像の推定の平均２乗誤差を、３つのデータ（ノイズ条件が無し，３０ｄｂ ＢＳＮＲ，２０ｄｂ ＢＳＮＲ）それぞれについて下表１及び図７に示す。本発明の超解像法を用いた推定結果の超解像度画像の方が、Tipping＆Bishopのガウス確率場を事前分布とした方法を用いた推定結果の超解像度画像よりも、精度が高く、平行移動量はいずれもサブピクセル精度での推定が実現できていることが観察される。【００７１】
【表１】
2007122838.tif 000021
【実施例２】
【００７２】
上述した実施例１では、事後分布の期待値（下式を参照）を高解像度画像の推定値としていた。【００７３】
【数１９】
2007122838.tif 000022
【００７４】
厳密な期待値計算には現実的な時間では実行不可能な積分計算がともなうため、計算コストをおさえるために適当な近似を行う必要がある。実施例１では近似に、事後分布を画素ごとに完全に独立なmixture of Gaussiansで置き換えるというものである。この近似により現実的な時間で推定値が計算可能になるが、近似による誤差のために、ノイズに引っ張られたざらついた高解像度画像を推定するという欠点が生じている。
実施例２では、事後分布の最大値（ＭＡＰ値）を推定画像と仮定し計算を行うものである。【００７５】
本実施例２では、ηについて和をとらず、ηもｘと同時に最適化を行っている。確率的な最適化手法であるシミュレーティッド・アニーリングを用いることで、モデルに近似を導入することなく、短時間での計算が可能となる。
図８が実施例２の実験に用いた原画像で、図９が１６枚の観測画像集合である。
図１０の（ａ）はノイズ無しの観測画像をもとに期待値推定を行った結果の高解像度画像で、（ｂ）はノイズ無しの観測画像をもとにＭＡＰ推定を行った結果の高解像度画像である。
図１１の（ａ）はＳＮ比が２０ｄＢの観測画像をもとに期待値推定を行った結果の高解像度画像で、（ｂ）は２０ｄＢの観測画像をもとにＭＡＰ推定を行った結果の高解像度画像である。
図１０と図１１を観察すると、ＭＡＰ推定は、期待値推定に比べて非常にはっきりした高解像度画像が得られていることが理解できる。
【産業上の利用可能性】
【００７６】
カメラ、ビデオなどは、その性能に応じて低解像度の画像を得るものであるが、そこから元の高解像度画像を復元することには大きなニーズがある。特に、同じ原画像から複数の劣化画像が得られることが仮定できる場合が多い。例えば手振れしながら撮像したビデオ画像の場合、時間局所的には同じ対象を撮ったものとみなすことができ、そうした連続した画像から元の対象を高精度で復元することができれば多大な応用になることは明らかである。他の例で言うならば、顕微鏡画像は、顕微鏡の物理的制約（例えば対物レンズのボケ過程、光の回折限界など）に応じて劣化が起こるが、多くの顕微鏡では、サンプリング周波数に応じて複数の撮影ができるため、それを用いて高精度で対象を復元できれば物理系の限界性能を超えた精度が実現でき得ることになる。本発明の超解像法および超解像プログラムは、こうした各種画像の復元に関するものであるため、多くの産業分野での応用が見込まれる。
【図面の簡単な説明】
【００７７】
【図１】本発明の超解像法のフローチャートを示した図
【図２】実施例１の実験に用いた原画像を示す。
【図３】実施例１の実験用に人工的に生成された１６枚の観測画像示す。
【図４】観測画像のうちの１枚を拡大した画像示す。
【図５】観測画像（図３に示したもの）に対し、本発明の超解像法を用いて推定を行ったノイズ無しのデータからの超解像度画像を（ａ）に、Tipping＆Bishopのガウス確率場を事前分布とした方法を用いて推定を行った超解像度画像を（ｂ）に示す。
【図６】観測画像（図３に示したもの）に対し、本発明の超解像法を用いて推定を行ったノイズ条件が３０ｄｂ ＢＳＮＲのデータからの超解像度画像を（ａ）に、Tipping＆Bishopのガウス確率場を事前分布とした方法を用いて推定を行った超解像度画像を（ｂ）に示す。
【図７】実施例１のデータ精度グラフ
【図８】実施例２の実験に用いた原画像を示す。
【図９】実施例２の実験用に１６枚の観測画像集合を示す。
【図１０】（ａ）はノイズ無しの観測画像をもとに期待値推定を行った結果の高解像度画像で、（ｂ）はノイズ無しの観測画像をもとにＭＡＰ推定を行った結果の高解像度画像を示している。
【図１１】（ａ）はＳＮ比が２０ｄＢの観測画像をもとに期待値推定を行った結果の高解像度画像で、（ｂ）は２０ｄＢの観測画像をもとにＭＡＰ推定を行った結果の高解像度画像を示している。
【図１】
2007122838.tif 000023
【図２】
2007122838.tif 000024
【図３】
2007122838.tif 000025
【図４】
2007122838.tif 000026
【図５】
2007122838.tif 000027
【図６】
2007122838.tif 000028
【図７】
2007122838.tif 000029
【図８】
2007122838.tif 000030
【図９】
2007122838.tif 000031
【図１０】
2007122838.tif 000032
【図１１】
2007122838.tif 000033
【国際調査報告】
2007122838.tif 000034
2007122838.tif 000035
2007122838.tif 000036
2007122838.tif 000037