JP2019003609A - 画像処理方法、画像処理装置、撮像装置および画像処理プログラム - Google Patents

Abstract

【課題】高速かつ高画質で被写体の画像を再構成可能な画像処理方法、画像処理装置、撮像装置および画像処理プログラムを提供すること。【解決手段】複数種類の符号化マスクを用いた被写体の撮像に基づく複数の符号化データを取得するステップと、複素数係数を用いた重み付き平均により、複数の符号化データに基づく合成データを算出するステップと、合成データおよび複数種類の符号化マスクに基づいて、被写体の画像を再構成処理によって算出するステップと、を有し、再構成処理は、複数種類の符号化マスクの透過率分布に基づいて、合成データに対し逆畳み込み演算を行うステップを有する。【選択図】図３

Description

本発明は、符号化マスクを用いた被写体を撮像することで生成された符号化データに基づく演算処理によって被写体の画像を再構成する画像処理方法に関する。

カメラや望遠鏡などの撮像装置は、レンズなどの光学素子を用いて被写体の像を撮像素子の撮像面に形成し、像の輝度分布を画像として取得する。撮像面上に形成された像は被写体の輝度分布をほぼそのまま再現しているため、撮像装置は特段の後処理なしに被写体の画像を取得することができる。しかしながら、理想的な画像を取得するためには高度に設計・加工されたレンズを使用する必要があり、そのようなレンズを備えた撮像装置は必然的に高価になってしまう。

特許文献１，２では、撮像装置の価格を抑えるために、レンズを用いることなく被写体の画像を再構成する方法が提案されている。上記方法では、撮像素子の前側（被写体側）に、特殊な形状の符号化マスクを配置し、取得された透過像から逆問題を計算機上で解くことで被写体の画像を再構成している。逆問題を解く方法として、特許文献１では、符号化マスクを、Ｔｏｅｐｌｉｔｚ行列を用いて分解し、ウィーナーフィルタおよびＬａｎｄｗｅｂｅｒ法を用いて再構成演算を行う方法が開示されている。また、特許文献２では、Ｔｉｋｈｏｎｏｖ正則化を用いた方法や圧縮検出を用いた方法などが開示されている。

米国特許第９４４５１１５号明細書 特表２０１６−５１０９１０号公報

通常、符号化マスクの周波数スペクトルには、振幅がゼロ近傍まで落ち込むゼロ点が存在する。そのため、周波数領域の情報を復元する場合、復元によってノイズが増幅されてしまい、良好な画像を算出することが困難である。一般的に、このようなノイズの増幅を抑えるために、ウィーナーフィルタを用いて回復処理が行われるが、回復処理はゼロ点を避けて行われるため、ゼロ点が存在する領域では実質的に回復効果が得られない。そのため、特許文献１の方法では、回復効果が得られない周波数成分が存在する。また、特許文献２の方法では、最適化問題を解くための計算負荷が大きく、画像を再構成するまでの時間が長くなってしまう。

本発明は、高速かつ高画質で被写体の画像を再構成可能な画像処理方法、画像処理装置、撮像装置および画像処理プログラムを提供することを目的とする。

本発明の一側面としての画像処理方法は、複数種類の符号化マスクを用いた被写体の撮像に基づく複数の符号化データを取得するステップと、複素数係数を用いた重み付き平均により、前記複数の符号化データに基づく合成データを算出するステップと、前記合成データおよび前記複数種類の符号化マスクに基づいて、前記被写体の画像を再構成処理によって算出するステップと、を有し、前記再構成処理は、前記複数種類の符号化マスクの透過率分布に基づいて、前記合成データに対し逆畳み込み演算を行うステップを有することを特徴とする。

また、本発明の他の側面としての画像処理装置は、複数種類の符号化マスクを用いた被写体の撮像に基づく複数の符号化データを取得する取得手段と、複素数係数を用いた重み付き平均により、前記複数の符号化データに基づく合成データを算出する算出手段と、前記合成データおよび前記複数種類の符号化マスクに基づいて、前記被写体の画像を再構成処理によって算出する再構成手段と、を有し、前記再構成手段は、前記再構成処理において、前記複数種類の符号化マスクの透過率分布に基づいて、前記合成データに対し逆畳み込み演算を行うことを特徴とする。

また、本発明の他の側面としての撮像装置は、撮像素子と、前記撮像素子に対して被写体側に配置された複数種類の符号化マスクと、前記撮像素子による前記複数種類の符号化マスクを用いた被写体の撮像に基づく複数の符号化データを取得する取得手段と、複素数係数を用いた重み付き平均により、前記複数の符号化データに基づく合成データを算出する算出手段と、前記合成データおよび前記複数種類の符号化マスクに基づいて、前記被写体の画像を再構成処理によって算出する再構成処理手段と、を有し、前記再構成処理手段は、前記再構成処理において、前記複数種類の符号化マスクの透過率分布に基づいて、前記合成データに対し逆畳み込み演算を行うことを特徴とする。

また、本発明の他の側面としての画像処理プログラムは、コンピュータに、複数種類の符号化マスクを用いた被写体の撮像に基づく複数の符号化データを取得するステップと、複素数係数を用いた重み付き平均により、前記複数の符号化データに基づく合成データを算出するステップと、前記合成データおよび前記複数種類の符号化マスクに基づいて、前記被写体の画像を再構成処理によって算出するステップと、を含む処理を実行させる画像処理プログラムであって、前記再構成処理は、複数種類の符号化マスクの透過率分布に基づいて、前記合成データに対し逆畳み込み演算を行うステップを含むことを特徴とする。

本発明によれば、高速かつ高画質で被写体の画像を再構成可能な画像処理方法、画像処理装置、撮像装置および画像処理プログラムを提供することができる。

本発明の実施形態に係る画像処理装置を有する撮像装置の概略図である。 被写体と符号化データとの関係図である。 本実施形態の画像処理方法を示すフローチャートである。 被写体を示す図である。 実施例１の符号化マスクを示す図である。 実施例１の符号化データを示す図である。 実施例１の再構成画像を示す図である。 従来方法を用いて取得された再構成画像を示す図である。 実施例２の符号化マスクを示す図である。 実施例２の再構成画像を示す図である。 従来方法を用いて取得された再構成画像を示す図である。 実施例３の符号化マスクを示す図である。 実施例３の再構成画像を示す図である。 従来方法を用いて取得された再構成画像を示す図である。 実施例４の符号化マスクを示す図である。 実施例４の再構成画像を示す図である。 従来方法を用いて取得された再構成画像を示す図である。 実施例６の画像処理方法を示すフローチャートである。 実施例６の再構成画像を示す図である。 符号化マスクの空間スペクトルの断面図である。

以下、本発明の好ましい実施の形態を、添付の図面に基づいて詳細に説明する。各図において、同一の部材については同一の参照番号を付し、重複する説明は省略する。

図１は、本発明の実施形態に係る画像処理装置１０４を有する撮像装置１００の概略図である。撮像装置１００は、符号化マスク１０１、撮像素子１０２、データ保持装置１０３および画像処理装置１０４を有する。

被写体１から発せられた光は、符号化マスク１０１を透過して撮像素子１０２に到達する。撮像素子１０２は、符号化マスク１０１によって符号化された光強度分布（符号化データ）を取得する。データ保持装置１０３は、撮像素子１０２により取得された符号化データを保存する。画像処理装置１０４は、データ保持装置１０３から読み込んだ符号化データに対して、マスク形状に対応する画像処理を行う。

以下、画像処理装置１０４が実行する画像処理方法について説明する。まず、図２を参照して、本発明の原理を説明する。図２は、被写体１と符号化データとの関係図である。

被写体１上の点から発せられた光は符号化マスク１０１をほぼ一様に照明し、撮像素子１０２の撮像面上には符号化マスク１０１の影が投影される。被写体１上の光軸Ｏ上の点でない点１ａから発せられた光は、符号化マスク１０１を斜めに照明する。そのため、撮像面上には、符号化マスク１０１に対して、光軸Ｏと物点１ａから発せられた光とがなす角θ、および符号化マスク１０１と撮像素子１０２との間隔Ｄによって定まるシフト量δだけずれた位置に符号化マスク１０１の影２０１が投影される。撮像素子１０２は、各物点から発せられた光によって形成された符号化マスク１０１の投影像が足し合わされたものを符号化データとして取得する。被写体１の輝度分布をずれ量δの関数ｏ（δ）、符号化マスク１０１の透過率を光軸Ｏに垂直な方向の位置座標ｘの関数ｔ（ｘ）とすると、撮像素子１０２によって取得される符号化データｉ（ｘ）は以下の式（１）で表すことができる。

ここで、＊は畳み込み積分を意味する。画像処理装置１０４は、式（１）を満足する被写体１の輝度分布ｏ（δ）を符号化マスク１０１の透過率ｔ（ｘ）および符号化データｉ（ｘ）から算出する。しかしながら、一般にこの演算は困難とされている。式（１）の両辺をフーリエ変換すると、畳み込みの定理により式（１）は以下の式（２）で表される。

Ｉ（ξ）＝Ｏ（ξ）×Ｔ（ξ） （２）
Ｉ（ξ），Ｏ（ξ），Ｔ（ξ）はそれぞれ符号化データｉ（ｘ）、被写体１の輝度分布ｏ（δ）、符号化マスク１０１の透過率ｔ（ｘ）のフーリエ変換であり、ξは空間周波数を表す。式（２）より、Ｉ（ξ）をＴ（ξ）で除算することでＯ（ξ）を取得することができる。しかしながら、一般に符号化マスクの周波数スペクトルＴは特定の周波数で値が非常に小さくなるゼロ点を有するため、除算によってその周波数領域での情報が過剰に増幅されたり、ノイズが増幅されてしまったりする。そのため、従来の方法では、ウィーナーフィルタ等を用いることで周波数領域でのゼロ落ちを避けていた。しかしながら、ウィーナーフィルタを用いた場合、回復効果が得られない周波数領域が生じるため、十分な画質の再生像を取得することができない。また、Ｌａｎｄｗｅｂｅｒ法などの繰り返し計算を用いる方法を使用する場合、演算に非常に重い負荷がかかってしまう。

本発明では、このような問題を回避するために、複数の符号化データを使用する。複数の符号化データを使用するためには、複数種類の符号化マスクを用意する必要がある。以下、一例として、第１の符号化マスクとしてｔ_１（ｘ）＝（１＋ｃｏｓ（２πξ_０ｘ^２））／２、第２の符号化マスクとしてｔ_２（ｘ）＝（１＋ｓｉｎ（２πξ_０ｘ^２））／２で表される透過率分布を有する２種類の符号化マスクを用意する。ξ_０は定数である。第１の符号化マスクおよび第２の符号化マスクを用いて取得される符号化データはそれぞれ、以下の式（３），（４）で表すことができる。

ｉ_１（ｘ）＝ｏ（ｘ）＊（１＋ｃｏｓ（２πξ_０ｘ^２））／２ （３）
ｉ_２（ｘ）＝ｏ（ｘ）＊（１＋ｓｉｎ（２πξ_０ｘ^２））／２ （４）
本発明では、画像処理装置１０４は、これらの符号化データに対して、以下の式（５）で表される重み付き平均を行うことで、合成データｉ_ｍｉｘ（ｘ）を算出する。

ｉ_ｍｉｘ（ｘ）＝ｉ_１＋ｊ×ｉ_２ （５）
ｊは虚数単位である。式（５）は、式（３）および式（４）を代入することで、以下の式（６）に書き換えることができる。

ｉ_ｍｉｘ（ｘ）＝ｏ（ｘ）＊ｅｘｐ（２πｊξ_０ｘ^２）／２＋Ｃ_１ （６）
Ｃ_１は、複素定数である。式（６）の第１項は、被写体の輝度分布と新たな複素関数ｅｘｐ（２πｊξ_０ｘ^２）との畳み込み演算である。式（６）は、複素定数Ｃ_１を除いて両辺をフーリエ変換すると、畳み込みの定理により以下の式（７）で表される。

Ｉ_ｍｉｘ（ξ）＝Ｏ（ξ）×Ｃ_２×ｅｘｐ（−πｊξ^２／２ξ_０） （７）
Ｃ_２は複素定数であり、Ｉ_ｍｉｘ（ξ）は合成データｉ_ｍｉｘ（ｘ）のフーリエ変換である。ｅｘｐ（−πｊξ^２／２ξ_０）は、全周波数領域で絶対値が１となる、すなわちゼロ点が存在しない関数である。したがって、合成データｉ_ｍｉｘ（ｘ）のスペクトルＩ_ｍｉｘ（ξ）をｅｘｐ（−πｊξ^２／２ξ_０）で除算する、すなわち周波数空間での逆フィルタを乗じることによって全ての周波数領域の情報を回復することができる。最後に、周波数領域の情報に対して、逆フーリエ変換すると、高画質な被写体の画像（再構成画像）を取得することができる。また、フーリエ変換および逆フーリエ変換は、ＦＦＴ（高速フーリエ変換）や逆ＦＦＴなどを用いることで、高速に実行することができる。

以上説明したように、上記方法を用いることで、複数の符号化データを式（５）に従って重み付き平均し、合成データを生成することで、逆フィルタによる再構成演算が可能となる。

上記説明では、特定の符号化マスクを使用したが、本発明はこれに限定されない。本発明の一般性を示すために、これまでの議論を一般化する。以下、任意の複素関数ｇ（ｘ）が任意の実関数Φ（ξ）を用いて以下の式（８），（９）で定義される場合について説明する。

Ｇ（ξ）＝ｅｘｐ（ｊΦ（ξ）） （８）
ｇ（ｘ）＝ＦＴ^−１［Ｇ（ξ）］ （９）
ここで、ＦＴはフーリエ変換を意味し、ＦＴ^−１は逆フーリエ変換を意味する。以下の説明では、Ｇ（ξ）を複素関数フィルタ、ｇ（ｘ）を複素関数マスクという。まず、複素関数マスクｇ（ｘ）の実部と虚部を用いて、以下の式（１０），（１１）で表される透過率分布を有する２種類の符号化マスクが作成される。

ｔ_１（ｘ）＝（１＋Ｒｅ［ｇ（ｘ）］）／２ （１０）
ｔ_２（ｘ）＝（１＋Ｉｍ［ｇ（ｘ）］）／２ （１１）
これらの符号化マスクを用いて取得される２つの符号化データ（符号化データ１および符号化データ２）に対して、式（５）と同様の複素数係数を使用した重み付き平均を行うことで、合成データｉ_ｍｉｘ（ｘ）は以下の式（１２）で表される。

ｉ_ｍｉｘ（ｘ）＝ｏ（ｘ）＊ｇ（ｘ）／２＋Ｃ_３ （１２）
Ｃ_３は、複素定数である。式（１２）の第１項は、式（６）と同様に、被写体の輝度分布と複素関数との畳み込み演算である。式（１２）は、複素定数Ｃ３を除いて両辺をフーリエ変換すると、畳み込みの定理により以下の式（１３）で表される。

Ｉ_ｍｉｘ（ξ）＝Ｏ（ξ）×Ｇ（ξ）／２ （１３）
式（８）よりＧ（ξ）は全周波数領域で絶対値が１となるため、Ｇ（ξ）で除算することによって全ての周波数領域の被写体のスペクトルを復元することができる。このように、本発明は、上述した一例のように実関数Φ（ξ）を−πξ^２／２ξ_０に特定した場合の符号化マスクの場合だけでなく、すなわち特定の符号化マスクに限定されることなく、高画質な被写体の画像を取得することができる。

なお、本実施形態では、簡単のため１次元の場合について説明したが、２次元の場合にも容易に拡張することができる。

以下、図３を参照して、画像処理装置１０４が実行する画像処理方法について説明する。図３は、本実施形態の画像処理方法を示すフローチャートである。本実施形態の画像処理方法は、ソフトウエアおよびハードウエア上で動作する画像処理プログラムにしたがって実行される。画像処理プログラムは、例えば、画像処理装置１０４内に格納されていてもよいし、コンピュータが読み取り可能な記録媒体に記録されていてもよい。また、本実施形態では画像処理装置１０４が画像処理方法を実行するが、パーソナルコンピュータ（ＰＣ）や専用の装置が画像処理装置として本実施形態の画像処理方法を実行してもよい。また、本実施形態の画像処理プログラムに対応する回路を設け、回路を動作させることで本実施形態の画像処理方法を実行してもよい。

ステップＳ１では、画像処理装置１０４は、複数種類の符号化マスクを用いた被写体１の撮像に基づく複数の符号化データを取得する。すなわち、画像処理装置１０４は、取得手段として機能する。

ステップＳ２では、画像処理装置１０４は、式（５）で表される、複素数係数を用いた重み付き平均により、ステップＳ１で取得した複数の符号化データに基づく合成データを算出する。すなわち、画像処理装置１０４は、算出手段として機能する。

ステップＳ３では、画像処理装置１０４は、ステップＳ２で算出した合成データからオフセットを減算する。

ステップＳ４では、画像処理装置１０４は、ＦＦＴを実行し、変換データ（スペクトルデータ）を算出する。

ステップＳ５では、画像処理装置１０４は、変換データに対して、複素関数フィルタＧ（ξ）から定まる周波数空間での逆フィルタを乗ずる。

ステップＳ６では、画像処理装置１０４は、逆ＦＦＴを実行することで、被写体１の画像を再構成する。

本実施形態では、ステップＳ３以降の処理を再構成処理といい、ステップＳ４からステップＳ６までの処理を逆畳み込み演算という。すなわち、画像処理装置１０４は、再構成手段として機能する。

以上の方法によって、複数の符号化データから被写体の画像を高速かつ高画質で再構成することができる。

以下の各実施例では、被写体として図４に示される画像を使用して再構成画像を取得する場合について説明する。

本実施例では、符号化マスク１０１として、図５に示される符号化マスクを使用する。本実施例の符号化マスクは、同心円状のパターンを有する。本実施例の符号化マスクは左の領域にｔ_１（ｘ）＝（１＋ｃｏｓ（２πξ_０（ｘ^２＋ｙ^２）））／２で表される透過率分布を有し、右の領域にｔ_２（ｘ）＝（１＋ｓｉｎ（２πξ_０（ｘ^２＋ｙ^２）））／２で表される透過率分布を有する。

まず、画像処理装置１０４は、図５に示される符号化マスクを用いて生成された図６に示される符号化データをデータ保持装置１０３から読み込む。図６に示される符号化データは、ξ_０＝０．０１７［１／ｐｉｘｅｌ］、サンプリング数として３００×６００を用いてシミュレーションによって生成される。

次に、画像処理装置１０４は、データ保持装置１０３から読み込んだ符号化データを左右２つの領域に分割し、各データをそれぞれ符号化データ１および符号化データ２として扱う。なお、画像処理装置１０４は、既に分割されてデータ保持装置１０３に保存された符号化データ１および符号化データ２を読み込んでもよい。

次に、画像処理装置１０４は、式（５）で表される重み付き平均により、符号化データ１および符号化データ２に基づく合成データを算出する。続いて、画像処理装置１０４は、合成データからオフセットを減算した後、ＦＦＴを実行することでスペクトルデータを算出する。

最後に、画像処理装置１０４は、スペクトルデータに対して、式（８）で定まる複素関数フィルタ（ｅｘｐ（−πｊ（ξ^２＋η^２）／２ξ_０））の逆数を乗じた後、逆ＦＦＴを実行する。ここで、ηは、ｙ方向の空間周波数である。

画像処理装置１０４は、本実施例の画像処理方法を実行することで、図７に示される、被写体を良く再現している被写体の画像を取得することができる。図８は、符号化データ１のみを用いて、本実施例と同等の計算時間でウィーナーフィルタにより再構成された被写体の画像を示している。図７と図８を比較すると、本実施例の画像処理方法を実行することで取得された画像がより高画質であることがわかる。画像評価指標の１つであるＰＳＮＲを算出すると、本実施例の画像処理方法を実行することで取得された画像では３４．７、従来方法を実行することで取得された画像では２６．８となり、定量的に見ても画質が向上している。

本実施例では、符号化マスク１０１として、図９に示される符号化マスクを使用する。本実施例では、実関数Φ（ξ，η）を−π（ξ^２＋０．５η^２）／２ξ_０とし、式（１０）を用いて算出された透過率分布を有する第１の符号化マスクおよび式（１１）を用いて算出された透過率分布を有する第２の符号化マスクが使用される。実関数Φは楕円形状を表すため、本実施例の符号化マスクは楕円形状を有する。

まず、画像処理装置１０４は、第１および第２の符号化マスクを用いて実施例１と同様のシミュレーションにより生成された、第１および第２の符号化マスクに対応する符号化データ１および符号化データ２を取得する。

次に、画像処理装置１０４は、式（５）で表される重み付き平均により、符号化データ１および符号化データ２に基づく合成データを算出する。続いて、画像処理装置１０４は、合成データからオフセットを減算した後、ＦＦＴを実行することでスペクトルデータを算出する。

最後に、画像処理装置１０４は、スペクトルデータに対して、式（８）で定まる複素関数フィルタの逆数を乗じた後、逆ＦＦＴを実行する。

画像処理装置１０４は、本実施例の画像処理方法を実行することで、図１０に示される、被写体を良く再現している被写体の画像を取得することができる。図１１は、符号化データ１のみを用いて、本実施例と同等の計算時間でウィーナーフィルタにより再構成された被写体の画像を示している。図１０と図１１を比較すると、本実施例の画像処理方法を実行することで取得された画像がより高画質であることがわかる。ＰＳＮＲを算出すると、本実施例の画像処理方法を実行することで取得された画像では２６．２、従来方法を実行することで取得された画像では２２．６となり、定量的に見ても画質が向上している。

本実施例では、符号化マスク１０１として、図１２に示される符号化マスクを使用する。本実施例では、実関数Φ（ξ，η）を−π（ξ^２＋０．２｜η｜）／（２ξ_０）とし、式（１０）を用いて算出された透過率分布を有する第１の符号化マスクおよび式（１１）を用いて算出された透過率分布を有する第２の符号化マスクが使用される。実関数Φは空間周波数の線形関数を有するため、本実施例の符号化マスクは概線状の形状を有する。

画像処理装置１０４は、図３のフローチャートに従って図１３に示される被写体の画像を取得することができる。図１４は、符号化データ１のみを用いて、本実施例と同等の計算時間でウィーナーフィルタにより再構成された被写体の画像を示している。図１３と図１４を比較すると、本実施例の画像処理方法を実行することで取得された画像がより高画質であることがわかる。ＰＳＮＲを算出すると、本実施例の画像処理方法を実行することで取得された画像では２１．８、従来方法を実行することで取得された画像では１７．１となり、定量的に見ても画質が向上している。

本実施例では、符号化マスク１０１として、図１５に示される符号化マスクを使用する。本実施例では、実関数Φ（ξ，η）を−π（｜ξ｜^３＋｜η｜^３）／ξ_０で表される３次関数とし、式（１０）で表される透過率分布を有する第１の符号化マスクおよび式（１１）で表される透過率分布を有する第２の符号化マスクが使用される。なお、本実施例では、実関数Φ（ξ，η）は、３次関数であるが、２次関数などのｎ次関数であってもよい。

画像処理装置１０４は、図３のフローチャートに従って図１６に示される被写体の画像を取得することができる。図１７は、符号化データ１のみを用いて、本実施例と同等の計算時間でウィーナーフィルタにより再構成された被写体の画像を示している。図１６と図１７を比較すると、本実施例の画像処理方法を実行することで取得された画像がより高画質であることがわかる。ＰＳＮＲを算出すると、本実施例の画像処理方法を実行することで取得された画像では２４．５、従来方法を実行することで取得された画像では２１．２となり、定量的に見ても画質が向上している。

本発明では繰り返し計算を伴わないため、高速な処理が可能であるが、演算の一部をあらかじめ解析的に実行しておくことで更に高速な処理を行うことも可能である。

本実施例では、実関数Φ（ξ）を−πξ^２／２ξ_０とする。簡単のため、１次元で結果を示す。周波数空間での逆フィルタが施された合成データはＦＴ^−１［Ｉ_ｍｉｘ（ξ）／（Ｃ_２×ｅｘｐ（−πｊξ^２／２ξ_０））］と表される。逆フィルタが施された合成データは、畳み込み定理により、以下の式（１４）で表される。

ｅｘｐ（４πｊξ_０ｘｘ’）を乗じて積分を行うことはフーリエ変換を行うことと等価であるため、この畳み込み演算は１回のＦＦＴで算出することができる。数値的に周波数空間での逆フィルタを乗ずるためには、少なくとも合成データのスペクトルを算出する場合と、フィルタが乗じられたスペクトルを逆フーリエ変換する場合の最低２回のＦＦＴを要する。しかしながら、本実施例で示したように、演算の一部をあらかじめ解析的に実行しておくことで、複素パターンを乗じた後、１回のＦＦＴで画像を再構成することが可能となり、演算の高速化を実現することができる。

合成データｉ_ｍｉｘ（ｘ）から被写体の輝度分布ｏ（ｘ）を算出する演算は、実空間での演算のみで実行することも可能である。以下の式（１５）は、式（１３）を変形したものである。

Ｏ（ξ）＝Ｉ_ｍｉｘ（ξ）／Ｇ（ξ）×２ （１５）
式（１５）は、両辺を逆フーリエ変換することで、以下の式（１６）で表される。

ｏ（ｘ）＝ｉ_ｍｉｘ（ｘ）＊ＦＴ^−１［２／Ｇ（ξ）］ （１６）
式（１６）において、ＦＴ^−１［２／Ｇ（ξ）］は実空間での逆フィルタである。式（１６）より、被写体の輝度分布は、合成データと実空間での逆フィルタとの畳み込み演算によって算出できることがわかる。実空間での逆フィルタは、被写体に依らないため、あらかじめ算出しておくことが可能である。そのため、画像を再構成する場合、合成データとの畳み込み演算のみを実行すればよい。

図１８は、画像処理装置１０４により実行される本実施例の画像処理方法を示すフローチャートである。

ステップＳ２１では、画像処理装置１０４は、複数種類の符号化マスクを用いた被写体１の撮像に基づく複数の符号化データを取得する。すなわち、画像処理装置１０４は、取得手段として機能する。

ステップＳ２２では、画像処理装置１０４は、式（５）で表される、複素数係数を用いた重み付き平均により、ステップＳ２１で取得した複数の符号化データに基づく合成データを算出する。すなわち、画像処理装置１０４は、算出手段として機能する。

ステップＳ２３では、画像処理装置１０４は、ステップＳ２２で算出した合成データからオフセットを減算する。

ステップＳ２４では、画像処理装置１０４は、ステップＳ２３でオフセットが減算された後の合成データとあらかじめ算出しておいた実空間での逆フィルタとの畳み込み演算（逆畳み込み演算）を行うことで、被写体１の画像を再構成する。

本実施例では、ステップＳ２４の処理を再構成処理という。すなわち、画像処理装置１０４は、再構成手段として機能する。

以上の方法によって、複数の符号化データから被写体の画像を高速かつ高画質で再構成することができる。

本実施例では、符号化マスク１０１として、実施例１と同じ図５に示される符号化マスクを使用する。画像処理装置１０４は、あらかじめこの符号化マスクに対応する実空間での逆フィルタを算出しておく。

画像処理装置１０４は、まず、第１および第２の符号化マスクを用いて実施例１と同様のシミュレーションにより生成された、第１および第２の符号化マスクに対応する符号化データ１および符号化データ２を取得する。次に、画像処理装置１０４は、式（５）で表される重み付き平均により、符号化データ１および符号化データ２に基づく合成データを算出する。続いて、画像処理装置１０４は、合成データからオフセットを減算する。最後に、画像処理装置１０４は、オフセットが除かれた合成データに対し、あらかじめ算出しておいた実空間での逆フィルタを畳み込む。図１９は、本実施例の画像処理方法を実行することで取得された被写体の画像であり、被写体を良く再現している。
［他の実施形態］
本発明は、上述の実施形態の１以上の機能を実現するプログラムを、ネットワーク又は記憶媒体を介してシステム又は装置に供給し、そのシステム又は装置のコンピュータにおける１つ以上のプロセッサーがプログラムを読出し実行する処理でも実現可能である。また、１以上の機能を実現する回路（例えば、ＡＳＩＣ）によっても実現可能である。

以上、本発明の好ましい実施形態について説明したが、本発明はこれらの実施形態に限定されず、その要旨の範囲内で種々の変形及び変更が可能である。

各実施例では、２つの符号化マスクを隣り合わせに配置する方法を説明したが、マスクの配置方法はこれに限られるものではない。空間的に異なる位置に異なる符号化マスクが配置されていればよく、更に言えば、複数の異なる符号化マスクによって符号化された符号化データを取得できれば本発明の効果を得ることができる。

本実施形態では、合成データから被写体の画像を再構成する方法として、合成データのスペクトルを複素関数フィルタＧ（ξ）で除算する方法を説明したが、これ以外の方法でも再構成は可能である。例えば、複素関数フィルタＧ（ξ）から生成されるウィーナーフィルタなどを合成データのスペクトルに乗じてもよい。逆畳み込み処理を行う演算方法であれば、本発明の効果を得ることができる。

本実施形態では、周波数スペクトルの値がゼロとなる箇所をゼロ点と称してきたが、ゼロ点として取り扱われるべき箇所は周波数スペクトルの値が厳密にゼロとなる箇所のみではない。特に、数値計算でデータを処理する場合、厳密なゼロ点というのは現れないことが多い。しかしながら、このような場合においても、他の領域に比べて周波数スペクトルの値が小さい領域は存在しており、その領域がゼロ点として称される箇所となる。例えば、ＦＦＴなどの数値計算で得られる周波数スペクトルの絶対値の最大値に対して、周波数スペクトルの値が十分に小さい点および領域はゼロ点として扱ってもよい。具体的には、周波数スペクトルの絶対値の最大値に対して、周波数スペクトルの値が０．１５倍以下の点および領域はゼロ点として扱ってもよい。また、式（１０）および式（１１）で定義される２つの符号化マスクの透過率には、負となる箇所を避けるために定数によるオフセットが足されている。このような符号化マスクの場合、周波数スペクトルの原点を除く領域での絶対値の最大値に対して、周波数スペクトルの値が０．１５倍以下の点および領域をゼロ点として扱えばよい。

本発明は、複数種類の符号化マスクのゼロ点が互いに重ならない場合に特に有効である。例えば、実関数Φ（ξ，η）をπ（ξ^２＋η^２）／０．０３４とした場合の第１および第２の符号化マスクの空間スペクトルの断面図を図２０に示す。実線が第１の符号化マスクの空間スペクトルの断面、破線が第２の符号化マスクの空間スペクトルの断面を表している。各符号化マスクの空間スペクトルにはゼロ点が存在するが、それぞれのゼロ点の位置はずれている。互いのゼロ点が重なっていないため、重み付き平均を行うことで、互いのマスクで取得できていない情報を補い合うことができ、結果として高画質で画像を再構成することができる。

本実施形態では、式（８）において、複素関数フィルタＧ（ξ）として全周波数領域で絶対値が１となる複素関数を用いているが、本発明はこれに限定されない。例えば、ａ（ξ、η）×ｅｘｐ（ｊΦ（ξ、η））で表される、絶対値に減衰や変調を有する関数であってもよい。特に、絶対値を表す関数ａ（ξ、η）がゼロ点を持たなければ、高い回復効果を得ることができる。また、厳密にゼロ点がなくなっている必要はない。ゼロ点が存在する場合は、ウィーナーフィルタ等を用いて回復処理を行えばよい。また、式（５）で表される重み付き平均により、ゼロ点の数またはゼロ点の領域の大きさが１つの符号化データを用いる場合に比べて小さくなっているので、１つの符号化データを用いた場合よりも高画質な再構成画像を取得することができる。

本実施形態では、空間周波数に対して特段の限定を行っていないが、実際に作成できるマスクには加工できる最小サイズが決まっており、これにより符号化マスクが有することが可能な最大周波数が決まってしまう。したがって、厳密に言えば、空間スペクトルにゼロ点が存在しない符号化マスクは存在しえない。しかしながら、本発明を実施する場合、対象とする周波数領域内で関数ａ（ξ、η）にゼロ点が存在しなければ、本発明の効果を得ることができる。例えば、以下の式（１７）で表される、符号化データのサンプリングピッチ（サンプリング間隔）Δの逆数から定まる周波数よりも小さな周波数領域内でゼロ点がなければよい。

ξ_ｍａｘ＝１／２Δ （１７）
本実施形態では、式（１０）および式（１１）を用いて、複素関数マスクｇ（ｘ）から符号化マスクを生成する方法を説明したが、本発明はこれに限定されない。例えば、複素関数マスクｇ（ｘ）の振幅が大きい場合、規格化などを施した後に式（１０）および式（１１）を用いて符号化マスクを生成してもよいし、式（１０）および式（１１）において加える定数を１とは異なる値で調整してもよい。各符号化マスクが複素関数ｇ（ｘ）の実部と虚部をそれぞれ含む形状をしていれば、本発明の効果を得ることができる。

更に言えば、第１および第２の符号化マスクの形状は、複素関数ｇ（ｘ）の実部と虚部を陽に含む形状である必要はない。実部と虚部から重み付き平均や線形変換等で生成される新たな関数１，２をそれぞれ含む符号化マスクを用いてもよい。この場合、合成データを生成する場合に、この重みから算出される複素数係数をかけて平均化すればよい。

本実施形態では、２種類の符号化マスクを用いる場合について説明したが、本発明はこれに限定されない。複素関数ｇ（ｘ）の実部と虚部を適当に分割して新たな関数を生成し、それらから３種類以上の符号化マスクを生成してもよい。

本実施形態では、実関数Φ（ξ）から生成された複素関数フィルタＧ（ξ）を用いて周波数空間および実空間での逆フィルタを算出したが、本発明はこれに限定されない。例えば、符号化マスクの透過率分布から複素数係数を用いた重み付き平均により、合成マスクｇ’（ｘ）を算出し、合成マスクｇ’（ｘ）からフーリエ変換によって合成スペクトルＧ’（ξ）を求め、合成スペクトルＧ’（ξ）から逆フィルタを生成してもよい。

本発明では、周波数スペクトルを算出する際、または周波数スペクトルから実空間でのデータを算出する際、フーリエ変換と逆フーリエ変換のどちらを用いてもよい。すべての演算において、一貫した組み合わせになっていればよい。変換方法の違いによって本発明の趣旨は変わるものではない。

式（１２）で算出される合成データが持つオフセットＣ_３を算出する方法は種々存在する。例えば、合成データの平均値をオフセットＣ_３として決定すればよい。また、符号化データの周辺領域の値をオフセットＣ_３として決定してもよい。また、画像を合成する前に符号化データからオフセットを除いても構わない。その際のオフセットとして、画像の平均値や、画像の周辺領域の値を用いればよい。

１０１ 符号化マスク

Claims (25)

1. 複数種類の符号化マスクを用いた被写体の撮像に基づく複数の符号化データを取得するステップと、
   複素数係数を用いた重み付き平均により、前記複数の符号化データに基づく合成データを算出するステップと、
   前記合成データおよび前記複数種類の符号化マスクに基づいて、前記被写体の画像を再構成処理によって算出するステップと、を有し、
   前記再構成処理は、前記複数種類の符号化マスクの透過率分布に基づいて、前記合成データに対し逆畳み込み演算を行うステップを有することを特徴とする画像処理方法。
2. 前記逆畳み込み演算では、前記複数種類の符号化マスクの透過率分布に基づいて算出される周波数空間での逆フィルタが用いられることを特徴とする請求項１に記載の画像処理方法。
3. 前記逆畳み込み演算は、
   フーリエ変換によって、前記合成データから変換データを算出するステップと、
   前記変換データに対して、前記周波数空間での逆フィルタを乗じるステップと、
   逆フーリエ変換を行うステップと、を有することを特徴とする請求項２に記載の画像処理方法。
4. 前記周波数空間での逆フィルタは、前記複数種類の符号化マスクの透過率分布から複素数係数を用いた重み付き平均により算出される合成マスクから、フーリエ変換によって算出される合成スペクトルに基づいて算出されることを特徴とする請求項２または３に記載の画像処理方法。
5. 前記逆畳み込み演算では、前記複数種類の符号化マスクの透過率分布に基づいて算出される実空間での逆フィルタが用いられることを特徴とする請求項１に記載の画像処理方法。
6. 前記逆畳み込み演算は、前記合成データに対して前記実空間での逆フィルタを畳み込むステップを有することを特徴とする請求項５に記載の画像処理方法。
7. 前記実空間での逆フィルタは、前記複数種類の符号化マスクの透過率分布から複素数係数を用いた重み付き平均により算出される合成マスクから、フーリエ変換によって算出される合成スペクトルに基づいて算出されることを特徴とする請求項５または６に記載の画像処理方法。
8. 前記再構成処理は、前記複数の符号化データから定数を減算するステップを有することを特徴とする請求項１から７のいずれか１項に記載の画像処理方法。
9. 前記再構成処理は、前記合成データから複素数の定数を減算するステップを有することを特徴とする請求項１から７のいずれか１項に記載の画像処理方法。
10. 前記定数は、前記複数の符号化データまたは前記合成データの平均値に基づいて決定されることを特徴とする請求項８または９に記載の画像処理方法。
11. 前記定数は、前記複数の符号化データまたは前記合成データの周辺領域の値に基づいて決定されることを特徴とする請求項８または９に記載の画像処理方法。
12. 前記複数種類の符号化マスクの透過率分布には、複数の関数の重み付き平均で表される第１の透過率分布と、前記第１の透過率分布とは異なる複数の関数の重み付き平均で表される第２の透過率分布と、が含まれることを特徴とする請求項１から１１のいずれか１項に記載の画像処理方法。
13. 前記第１の透過率分布の周波数スペクトルの最大値に対して０．１５倍以下となる領域と、前記第２の透過率分布の周波数スペクトルの最大値に対して０．１５倍以下となる領域と、は異なることを特徴とする請求項１２に記載の画像処理方法。
14. 前記複素数係数は、前記第１の透過率分布を表す複数の関数の重み付き平均の重みおよび前記第２の透過率分布を表す複数の関数の重み付き平均の重みに基づいて算出されることを特徴とすることを特徴とする請求項１２または１３に記載の画像処理方法。
15. 前記複数種類の符号化マスクの透過率分布は、複素関数で表される複素関数マスクに基づいて生成されることを特徴とする請求項１２から１４のいずれか１項に記載の画像処理方法。
16. 前記第１および第２の透過率分布は、前記複素関数に基づくことを特徴とする請求項１５に記載の画像処理方法。
17. 前記第１および第２の透過率分布は、前記複素関数の実部および虚部の重み付き平均により表され、
    前記複素関数は、複素関数フィルタの逆フーリエ変換で表され、
    前記逆畳み込み演算は、
    フーリエ変換によって、前記合成データから変換データを算出するステップと、
    前記複素関数フィルタに基づいて算出される周波数空間での逆フィルタを乗じるステップと、
    逆フーリエ変換を行うステップと、を有することを特徴とする請求項１５または１６に記載の画像処理方法。
18. 前記第１および第２の透過率分布は、前記複素関数の実部および虚部の重み付き平均により表され、
    前記複素関数は、複素関数フィルタの逆フーリエ変換で表され、
    前記逆畳み込み演算は、前記複素関数フィルタに基づいて算出される実空間での逆フィルタを前記合成データに畳み込むステップを有することを特徴とする請求項１５または１６に記載の画像処理方法。
19. 前記周波数空間での複素関数フィルタは、２次元の各方向の空間周波数をそれぞれξ，ηとするとき、正の実関数ａ（ξ,η）および実関数Φ（ξ,η）を用いて、ａ（ξ,η）×ｅｘｐ（ｊΦ（ξ,η））と表され、
    前記正の実関数ａ（ξ,η）は、前記符号化データのサンプリング間隔をΔ、前記正の実関数ａ（ξ,η）の最大値をａ_ｍａｘとするとき、
    （ξ^２＋η^２）^１／２≦１／２Δ
    なる条件式を満足する領域内で、
    ａ（ξ,η）≧０．１５ａ_ｍａｘ
    なる条件式を満足することを特徴とする請求項１７に記載の画像処理方法。
20. 前記実関数Φ（ξ,η）は、前記空間周波数ξ，ηの少なくとも一方の線形関数を含むことを特徴とする請求項１９に記載の画像処理方法。
21. 前記実関数Φ（ξ,η）は、前記空間周波数ξ，ηの少なくとも一方の２次関数を含むことを特徴とする請求項１９に記載の画像処理方法。
22. 前記実関数Φ（ξ,η）は、前記空間周波数ξ，ηの少なくとも一方の３次関数を含むことを特徴とする請求項１９に記載の画像処理方法。
23. 複数種類の符号化マスクを用いた被写体の撮像に基づく複数の符号化データを取得する取得手段と、
    複素数係数を用いた重み付き平均により、前記複数の符号化データに基づく合成データを算出する算出手段と、
    前記合成データおよび前記複数種類の符号化マスクに基づいて、前記被写体の画像を再構成処理によって算出する再構成手段と、を有し、
    前記再構成手段は、前記再構成処理において、前記複数種類の符号化マスクの透過率分布に基づいて、前記合成データに対し逆畳み込み演算を行うことを特徴とする画像処理装置。
24. 撮像素子と
    前記撮像素子に対して被写体側に配置された複数種類の符号化マスクと、
    前記撮像素子による前記複数種類の符号化マスクを用いた被写体の撮像に基づく複数の符号化データを取得する取得手段と、
    複素数係数を用いた重み付き平均により、前記複数の符号化データに基づく合成データを算出する算出手段と、
    前記合成データおよび前記複数種類の符号化マスクに基づいて、前記被写体の画像を再構成処理によって算出する再構成処理手段と、を有し、
    前記再構成処理手段は、前記再構成処理において、前記複数種類の符号化マスクの透過率分布に基づいて、前記合成データに対し逆畳み込み演算を行うことを特徴とする撮像装置。
25. コンピュータに、
    複数種類の符号化マスクを用いた被写体の撮像に基づく複数の符号化データを取得するステップと、
    複素数係数を用いた重み付き平均により、前記複数の符号化データに基づく合成データを算出するステップと、
    前記合成データおよび前記複数種類の符号化マスクに基づいて、前記被写体の画像を再構成処理によって算出するステップと、を含む処理を実行させる画像処理プログラムであって、
    前記再構成処理は、複数種類の符号化マスクの透過率分布に基づいて、前記合成データに対し逆畳み込み演算を行うステップを含むことを特徴とする画像処理プログラム。