(57)【要約】【課題】高解像度化処理により増幅されたノイズを低減する。【解決手段】画像処理方法は、光学系を用いた撮像により力画像を取得し（スＳ１０１）、該入力画像に含まれる第１のノイズ量を取得し（Ｓ１０２）、入力画像に対して高解像度化処理を行うことにより第１の画像を生成し（Ｓ１０５）、該第１の画像をそれぞれ異なる周波数帯域成分としての複数の第２の画像に周波数分解し（Ｓ１０６）、第１のノイズ量と高解像度化処理の周波数特性とに基づいて、各第２の画像において複数に分割された領域ごとに第２のノイズ量を算出し（Ｓ１０８）、該第２のノイズ量に基づいて、各第２の画像の該領域ごとにノイズ低減処理を行い（Ｓ１０９）、該ノイズ低減処理が行われた複数の第２の画像を再合成して第３の画像を生成する（Ｓ１１０）。【選択図】図１
2015166979.tif 000002
【特許請求の範囲】
【請求項１】
光学系を用いた撮像により入力画像を取得し、
該入力画像に含まれる第１のノイズ量を取得し、
前記入力画像に対して高解像度化処理を行うことにより第１の画像を生成し、
該第１の画像をそれぞれ異なる周波数帯域成分としての複数の第２の画像に周波数分解し、
前記第１のノイズ量と前記高解像度化処理の周波数特性とに基づいて、前記各第２の画像において複数に分割された領域ごとに第２のノイズ量を算出し、
該第２のノイズ量に基づいて、前記各第２の画像の前記領域ごとにノイズ低減処理を行い、
該ノイズ低減処理が行われた複数の第２の画像を再合成して第３の画像を生成することを特徴とする画像処理方法。
【請求項２】
前記高解像度化処理は、前記光学系の光学特性に基づいて行われることを特徴とする請求項１に記載の画像処理方法。
【請求項３】
前記ノイズ低減処理は、前記複数の第２の画像において前記光学系の光学特性に基づいて分割された前記領域ごとに行われることを特徴とする請求項１または２に記載の画像処理方法。
【請求項４】
前記第１のノイズ量は、前記複数の第２の画像のうち最も高周波側の第２の画像において前記光学系の光学特性に基づいて分割された前記領域ごとに取得されることを特徴とする請求項１から３のいずれか一項に記載の画像処理方法。
【請求項５】
前記第１のノイズ量は、前記撮像の条件に対するノイズ量を測定して得られたノイズ量データに基づいて、前記光学系の光学特性に基づいて分割された前記領域ごとの平均輝度値を用いて決定することを特徴とする請求項１から３のいずれか一項に記載の画像処理方法。
【請求項６】
前記撮像の条件は、ＩＳＯ感度と輝度レベルであることを特徴とする請求項５に記載の画像処理方法。
【請求項７】
前記第２のノイズ量は、前記光学系の光学特性に基づいて分割された前記領域ごとの前記高解像度化処理の周波数特性と前記周波数分解に用いる基底関数の周波数特性とに基づいて算出されることを特徴とする請求項１から６のいずれか一項に記載の画像処理方法。
【請求項８】
前記第２のノイズ量は、前記光学系の光学特性に基づいて分割された前記領域ごとの前記高解像度化処理の周波数特性と前記周波数分解に用いる基底関数の周波数特性とに基づいて算出された、前記高解像度化処理における前記複数の領域において一様ではない増幅量に基づいて算出されることを特徴とする請求項１から６のいずれか一項に記載の画像処理方法。
【請求項９】
前記周波数分解は、離散ウェーブレット変換により行われることを特徴とする請求項１，４，７および８のいずれか一項に記載の画像処理方法。
【請求項１０】
光学系を用いた撮像により生成された入力画像を取得する手段と、
該入力画像に含まれる第１のノイズ量を取得する手段と、
前記入力画像に対して高解像度化処理を行うことにより第１の画像を生成する手段と、
該第１の画像をそれぞれ異なる周波数帯域成分としての複数の第２の画像に周波数分解する手段と、
前記第１のノイズ量と前記高解像度化処理の周波数特性とに基づいて、前記各第２の画像において複数に分割された領域ごとに第２のノイズ量を算出する手段と、
該第２のノイズ量に基づいて、前記各第２の画像の前記領域ごとにノイズ低減処理を行う手段と、
該ノイズ低減処理が行われた複数の第２の画像を再合成して第３の画像を生成する手段とを有することを特徴とする画像処理装置。
【請求項１１】
光学系を用いた撮像を行って入力画像を生成する手段と、
請求項１０に記載の画像処理装置とを有することを特徴とする撮像装置。
【請求項１２】
画像処理が可能なコンピュータに、
光学系を用いた撮像により入力画像を取得し、
該入力画像に含まれる第１のノイズ量を取得し、
前記入力画像に対して高解像度化処理を行うことにより第１の画像を生成し、
該第１の画像をそれぞれ異なる周波数帯域成分としての複数の第２の画像に周波数分解し、
前記第１のノイズ量と前記高解像度化処理の周波数特性とに基づいて、前記各第２の画像において複数に分割された領域ごとに第２のノイズ量を算出し、
該第２のノイズ量に基づいて、前記各第２の画像の前記領域ごとにノイズ低減処理を行い、
該ノイズ低減処理が行われた複数の第２の画像を再合成して第３の画像を生成する前記画像処理を実行させるコンピュータプログラムであることを特徴とする画像処理プログラム。
【発明の詳細な説明】
【技術分野】
【０００１】
本発明は、撮像により生成された画像のノイズを低減する画像処理技術に関し、特に高解像度化処理を行った画像のノイズを低減する技術に関する。
【背景技術】
【０００２】
被写体の１点から出た光は、光学系を通ると、該光学系の回折や収差の影響によって１点に収束することができず、微小な広がりを持つ。このような微小な広がりを持った分布を点像分布関数（Point Spread Function：ＰＳＦ）と呼ぶ。撮像光学系を用いた撮像により生成された画像は、理想的な被写体像にこのＰＳＦが畳み込まれたものであるため、該画像の解像度が劣化する（ぼける）。
【０００３】
このため、このような劣化画像を補正する高解像度化処理が行われる場合がある。高解像度化処理の１つとして、画像回復処理がある。画像回復処理は、例えば、撮像光学系の光学伝達関数（Optical Transfer Function：ＯＴＦ）の逆特性を有する画像回復フィルタを、補正対象画像（入力画像）に対して畳み込むこと（コンボリューション）により行われる。ＯＴＦはＰＳＦをフーリエ変換することで得られる。
【０００４】
ただし、画像回復処理を行うことで、撮像に用いられた撮像装置（撮像素子）に起因して劣化画像に含まれていたノイズを増幅させたりリンギングを発生させたりするおそれがある。このため、単純にＯＴＦの逆特性を有する画像回復フィルタを用いた画像回復処理による高解像度化処理では、良好な高解像度画像を得ることは困難である。このような問題を解消することが可能な画像回復処理方法として、画像信号とノイズ信号の強度比（ＳＮＲ）に応じて回復度合いを制御可能なウィナーフィルタを用いる方法が知られている。ウィナーフィルタには、ノイズの増幅率を制御する調整項（増幅率調整項）が設けられている。
【０００５】
ところが、ウィナーフィルタにおいて増幅率調整項を増幅率が強くなるように設定すると、画像回復による高解像度化の効果は向上するが、ノイズが著しく増幅されてしまう。一方、増幅率調整項を増幅率が弱くなるように設定すると、ノイズの増幅は抑制できるが、高解像度化の効果も低下する。
【０００６】
そこで、高解像度化処理により生成された高解像度画像に対してノイズ低減処理を行うことで高解像度のノイズが少ない画像を取得する方法として、ウェーブレット変換を用いる方法が特許文献１や非特許文献１にて開示されている。特許文献１には、高解像度化処理として振れ補正を行った場合において、ウェーブレット空間上でノイズを低減する方法が開示されている。また、非特許文献１には、画像全体に一様に高解像度化処理としての画像回復処理を行った場合において、ウェーブレット空間上にて、高解像度化処理の周波数特性に基づいて画像全体に一様にノイズ低減処理を行い、増幅されたノイズを低減する方法が開示されている。
【先行技術文献】
【特許文献】
【０００７】
【特許文献１】米国特許公開２００７−０１６５９６１号公報
【非特許文献】
【０００８】
【非特許文献１】R. Neelamani et al.”ForWaRD: Fourier-Wavelet Regularized Deconvolution for Ill-Conditioned Systems” IEEE Transactions on signal processing, vol.52, no.2 (2004).
【発明の概要】
【発明が解決しようとする課題】
【０００９】
しかしながら、特許文献１での高解像度化処理としての振れ補正は、画像全体に一様に生じた振れを補正対象としているのに対して、撮像光学系による劣化は画像全体において一様ではなく、それに対する高解像度化処理は画像内の領域ごとに行う必要がある。特許文献１には、画像内の領域ごとに高解像度化処理を行った場合において増幅されたノイズを低減する方法については開示されていない。
【００１０】
また、画像内の領域ごとに高解像度化処理を行った場合において、非特許文献１にて開示されているように画像全体に一様にノイズ低減処理を行うと、領域ごとの高解像度化処理によるノイズの増幅量に応じて残留ノイズやぼけが発生する。
【００１１】
本発明は、画像の領域ごとに高解像度化処理を行った場合であっても、高解像度化処理により増幅されたノイズを効果的に低減することができるようにした画像処理方法、画像処理装置、撮像装置および画像処理プログラムを提供する。
【課題を解決するための手段】
【００１２】
本発明の一側面としての画像処理方法は、光学系を用いた撮像により入力画像を取得し、該入力画像に含まれる第１のノイズ量を取得し、入力画像に対して高解像度化処理を行うことにより第１の画像を生成し、該第１の画像をそれぞれ異なる周波数帯域成分としての複数の第２の画像に周波数分解し、第１のノイズ量と高解像度化処理の周波数特性とに基づいて、各第２の画像において複数に分割された領域ごとに第２のノイズ量を算出し、該第２のノイズ量に基づいて、各第２の画像の該領域ごとにノイズ低減処理を行い、該ノイズ低減処理が行われた複数の第２の画像を再合成して第３の画像を生成することを特徴とする。
【００１３】
また、本発明の他の一側面としての画像処理装置は、光学系を用いた撮像により入力画像を取得する手段と、該入力画像に含まれる第１のノイズ量を取得する手段と、入力画像に対して高解像度化処理を行うことにより第１の画像を生成する手段と、該第１の画像をそれぞれ異なる周波数帯域成分としての複数の第２の画像に周波数分解する手段と、第１のノイズ量と高解像度化処理の周波数特性とに基づいて、各第２の画像において複数に分割された領域ごとに第２のノイズ量を算出する手段と、該第２のノイズ量に基づいて、各第２の画像の該領域ごとにノイズ低減処理を行う手段と、該ノイズ低減処理が行われた複数の第２の画像を再合成して第３の画像を生成する手段とを有することを特徴とする。
【００１４】
なお、光学系を用いた撮像を行って入力画像を生成する手段と、上記画像処理装置とを有する撮像装置も、本発明の他の一側面を構成する。
【００１５】
さらに、本発明の他の一側面としての画像処理プログラムは、画像処理が可能なコンピュータに、光学系を用いた撮像により生成された入力画像を取得し、該入力画像に含まれる第１のノイズ量を取得し、入力画像に対して高解像度化処理を行うことにより第１の画像を生成し、該第１の画像をそれぞれ異なる周波数帯域成分としての複数の第２の画像に周波数分解し、第１のノイズ量と高解像度化処理の周波数特性とに基づいて、各第２の画像において複数に分割された領域ごとに第２のノイズ量を算出し、該第２のノイズ量に基づいて、各第２の画像の該領域ごとにノイズ低減処理を行い、該ノイズ低減処理が行われた複数の第２の画像を再合成して第３の画像を生成する前記画像処理を実行させるコンピュータプログラムであることを特徴とする。
【発明の効果】
【００１６】
本発明では、高解像度化処理により生成された第１の画像を周波数分解して得られた複数の第２の画像における各第２の画像の領域ごとにノイズ低減処理を行い、ノイズ低減処理が行われた複数の第２の画像を再合成して第３の画像を生成する。このため、本発明によれば、入力画像に対して領域ごとの高解像度化処理が行われた場合でも、高解像度化処理により増幅されたノイズを効果的に低減することができる。
【図面の簡単な説明】
【００１７】
【図１】本発明の実施例１である画像処理方法を示すフローチャート。
【図２】実施例１における画像回復処理のシフトインバリアント領域と光学伝達関数を示す図。
【図３】実施例１における高解像度化処理を行った画像をウェーブレット変換により複数のサブバンド画像に周波数分解した一例を示す図。
【図４】実施例１におけるサブバンド画像の領域を示す図。
【図５】実施例１におけるサブバンド画像の領域ごとの増幅量を示す図。
【図６】本発明の実施例２において行うノイズ低減処理を示すフローチャート。
【図７】本発明の実施例３である画像処理装置を内蔵した撮像装置を示す図。
【発明を実施するための形態】
【００１８】
以下、本発明の実施例について図面を参照しながら説明する。
【００１９】
まず、本発明の前提技術について説明する。撮像光学系（以下、単に光学系ともいう）を用いた撮像により得られた撮影画像（以下、入力画像ともいう）における該光学系の回折や収差に起因する劣化を補正するための高解像度化処理の１つとして、以下の画像回復処理が行われる。
【００２０】
実空間（ｘ，ｙ）上において光学系による劣化を受ける前の画像をｆ（ｘ，ｙ）とし、光学系のＰＳＦをｈ（ｘ，ｙ）とし、劣化した画像をｇ（ｘ，ｙ）とすると、これらの関係は以下に示す（１）式で表現することができる。
ｇ（ｘ，ｙ）＝∫∫ｆ（Ｘ，Ｙ）・ｈ（ｘ−Ｘ，ｙ−Ｙ）ｄＸｄＹ  ・・・（１）
そして、（１）式に対してフーリエ変換を行い、実空間（ｘ，ｙ）から周波数空間（ｕ，ｖ）への変換を行うと、
Ｇ（ｕ，ｖ）＝Ｆ（ｕ，ｖ）・Ｈ（ｕ，ｖ）  ・・・（２）
となる。Ｆ（ｕ，ｖ）はｆ（ｘ，ｙ）のフーリエ変換であり、Ｇ（ｕ，ｖ）はｇ（ｘ，ｙ）のフーリエ変換であり、Ｈ（ｕ，ｖ）はｈ（ｘ，ｙ）のフーリエ変換である。すなわち、Ｈ（ｕ，ｖ）はＯＴＦ（光学伝達関数）である。
【００２１】
さらに、（２）式を変形すると、
Ｆ（ｕ，ｖ）＝Ｇ（ｕ，ｖ）／Ｈ（ｕ，ｖ）  ・・・（３）
となる。これは周波数空間上で劣化画像ｇ（ｘ，ｙ）のフーリエ変換Ｇ（ｕ，ｖ）をＰＳＦであるｈ（ｘ，ｙ）のフーリエ変換Ｈ（ｕ，ｖ）で除することにより、劣化を受ける前の画像ｆ（ｘ，ｙ）のフーリエ変換であるＦ（ｕ，ｖ）を得ることができることを意味する。したがって、Ｆ（ｕ，ｖ）に逆フーリエ変換を行うことにより、劣化前の画像ｆ（ｘ，ｙ）を得ることができる。
【００２２】
しかし、実際にこのような処理を行って劣化前の画像を得ようとすると、撮像素子によって生じるノイズを著しく増幅させたり画像にリンギングを発生させたりする可能性が高い。このため、このような画像処理による弊害の発生により、単純な逆特性で除算する手法では良好な画像を得ることが困難である。そこで、ノイズの増幅を抑制するための画像回復方法として、（４）式で表されるウィナーフィルタを用いることが知られている。
【００２３】
上述したようにＨ（ｕ，ｖ）は光学系のＯＴＦであり、Γはノイズの増幅量を低減するための定数である。換言すれば、画像回復処理の対象画像に全くノイズが発生していなくても、増幅率が大きくなりすぎると、リンギングやその他のエッジの掘り込み（アンダーシュート）のような弊害が発生し易くなる。このため、Γは増幅率を制御する調整項として利用することが可能であると言える。また、Γを定数とせずに、Γ（ｕ，ｖ）として、周波数特性を考慮した調整項としてもよい。以下の説明では、λｆ（ｕ，ｖ）を高解像度化処理の周波数特性といい、Γを増幅率調整項という。
【００２４】
【数１】
2015166979.tif 000003
【００２５】
また、光学系の光学伝達関数Ｈ（ｕ，ｖ）は像面、つまりは入力画像における領域（以下、画像領域という）ごとに異なる。したがって、（４）式で示す高解像度化処理の周波数特性λf（ｕ，ｖ）も画像領域ごとに異なる特性としてもよい。ただし、ここでいう画像領域とは、画素単位の領域も含む。
【００２６】
（４）式を、劣化画像ｇ（ｘ，ｙ）のフーリエ変換Ｇ（ｕ，ｖ）に乗算すれば、光学系の回折や収差によって発生したＰＳＦの位相成分を０とし、振幅成分の周波数特性を増幅することで高解像度かつ良好な画像を得ることができる。すなわち、ウィナーフィルタを用いた画像回復処理により回復された画像の周波数空間情報をＲ（ｕ，ｖ）とすると、
Ｒ（ｕ，ｖ）＝Ｇ（ｕ，ｖ）・λｆ（ｕ，ｖ）  ・・・（５）
が得られる。
【００２７】
ウィナーフィルタにおいて、Γを調節して増幅率を強くすると、高解像度化処理の効果は向上するが、ノイズが著しく増幅されてしまう。また、Γを調節して増幅率を弱くすると、ノイズの増幅を抑制することができるが、同時に高解像度化処理の効果が低下してしまう。
【００２８】
そこで、高解像度化処理を行って取得した高解像度画像に対してノイズ低減処理を行うことによって、高解像度かつ良好な画像を取得できる。その方法として、前述したウェーブレット変換を用いた方法が知られている。
【００２９】
ウェーブレット変換とは、局在する小さな波（ウェーブレット）を用いて画像の位置ごとに周波数解析を行い、信号を高周波成分と低周波成分に分解していく手法である。画像位置ごとに周波数解析を行うため、周波数領域と空間領域の両方で情報の表現が可能であり、画像位置ごとに異なる処理に対して有効である。
【００３０】
ウェーブレット変換では、（６）式のように、信号ｆ（ｘ）と空間的に局在する基底関数ψ（ｘ）との畳み込み積分を行う。
【００３１】
【数２】
2015166979.tif 000004
【００３２】
（７）式に示す基底関数ψ（ｘ）には、スケールパラメータａとシフトパラメータｂとが含まれており、それぞれ周波数方向および空間方向の位置に対応している。
【００３３】
【数３】
2015166979.tif 000005
【００３４】
離散ウェーブレット変換は、（８）式のようなスケールパラメータａとシフトパラメータｂとを離散化した基底関数ψ（ｘ）を用いることで行われる。
【００３５】
【数４】
2015166979.tif 000006
【００３６】
ただし、ｊはウェーブレット変換レベルであり、ｋは位置パラメータを示す。
【００３７】
また、（９）式に示す離散ウェーブレット変換では、正規直交基底を構成するスケーリング関数φとウェーブレット関数ψとの組からなる基底関数を用いて情報を表現する。スケーリング関数φとウェーブレット関数ψは数列で与えられ、スケーリング数列ｐｋとウェーブレット数列ｑｋを用いて情報を表現する。
【００３８】
【数５】
2015166979.tif 000007
【００３９】
ただし、ｓｊ，ｋは低周波成分であるスケーリング係数であり、ｗｊ，ｋは高周波成分であるウェーブレット係数である。基底関数については多数提案されており、使用したい関数を選択し、そのスケーリング数列ｐｋを用いてウェーブレット変換を行う。ウェーブレット数列ｑｋは、スケーリング数列ｐｋを用いて、（１０）式のようなツースケール関係によって導くことができる。
【００４０】
【数６】
2015166979.tif 000008
【００４１】
画像のウェーブレット変換では、画像の水平方向に対してウェーブレット変換を行って低周波成分と高周波成分に分解し、さらに分解して得られた低周波成分と高周波成分の垂直方向に対してウェーブレット変換を行う。ウェーブレット変換により、画像は４分割され、周波数帯域が互いに異なる４つのサブバンド画像に周波数分解される。このとき、左上の低周波帯域成分（スケーリング係数）のサブバンド画像をＬＬ１とし、右下の高周波帯域成分（ウェーブレット係数）のサブバンド画像をＨＨ１という。また、右上（ＨＬ１）と左下（ＬＨ１）のサブバンド画像はそれぞれ、水平方向に高周波帯域成分をとって垂直方向に低周波帯域成分を取り出したもの、および水平方向に低周波帯域成分をとって垂直方向に高周波帯域成分を取り出したものである。
【００４２】
さらに、ＬＬ１をウェーブレット変換すると、画像サイズを半分にして、ＬＬ２，ＨＬ２，ＬＨ２，ＨＨ２のサブバンド画像に分解することができ、分解して得られたＬＬに対して変換レベル回数だけ分解することができる。
【００４３】
ウェーブレット変換を用いてノイズ低減処理を行う方法として、しきい値処理が知られている。しきい値処理とは、設定したしきい値より小さい量の成分をノイズとみなし、そのノイズを低減する方法である。ウェーブレット空間上でのしきい値処理は、ＬＬ以外のサブバンド画像に対して行われ、（１１）式に示すように、しきい値以下の絶対値を持つウェーブレット係数ｗｓｕｂｂａｎｄ（ｘ，ｙ）を０に置き換える。これにより、ノイズ成分が低減される。
【００４４】
【数７】
2015166979.tif 000009
【００４５】
上記しきい値は高解像度画像に含まれるノイズ量に基づいて設定され、具体的には標準偏差が用いられる。ρｓｕｂｂａｎｄは重みパラメータである。
高解像度画像をウェーブレット変換により周波数分解することで得られたサブバンド画像におけるノイズ量σｓｕｂｂａｎｄは、（１２）式から求められる。すなわち、高解像度化処理を行う前の入力画像に含まれるノイズ量σ０と、高解像度化処理によるサブバンド画像ごとの増幅量Ａｓｕｂｂａｎｄとから求めることができる。
【００４６】
【数８】
2015166979.tif 000010
【００４７】
【数９】
2015166979.tif 000011
【００４８】
入力画像に含まれるノイズ量σ０は、該入力画像から測定または推定することにより求められる。ノイズが実空間および周波数空間において一様なホワイトガウシアンノイズであるとき、入力画像に含まれるノイズを、（１４）式に示すようなＭＡＤ（Median Absolute Deviation）から推定する手法が知られている。
ＭＡＤ＝ｍｅｄｉａｎ（｜ｗＨＨ１−ｍｅｄｉａｎ（ｗＨＨ１）｜）    ・・・（１４）
ＭＡＤは入力画像をウェーブレット変換して得られたＨＨ１のサブバンド画像におけるウェーブレット係数ｗＨＨ１のメディアン（中央値）を用いて求められ、標準偏差とＭＡＤが（１５）式の関係であることからノイズ成分の標準偏差を推定することができる。
σ０＝ＭＡＤ／０．６７４５  ・・・（１５）
高解像度化処理によるサブバンド画像ごとの増幅量Ａｓｕｂｂａｎｄは、（１３）式より求められ、高解像度化処理の周波数特性λｆ（ｕ，ｖ）に基づいて決定される。上記したように、高解像度化処理の周波数特性λｆ（ｕ，ｖ）は、光学伝達関数Ｈ（ｕ，ｖ）から決定される。Ψｓｕｂｂａｎｄ（ｕ，ｖ）はウェーブレット変換による周波数分解の周波数特性であり、ウェーブレット変換時に使用した基底関数から一意に決定される。Ｎはλｆ（ｕ，ｖ）のデータ数である。
【００４９】
上記の方法により算出した高解像度画像をウェーブレット変換により分解したサブバンド画像のノイズ量σｓｕｂｂａｎｄは、高解像度化処理による増幅量Ａｓｕｂｂａｎｄに基づいたものである。つまり、該ノイズ量σｓｕｂｂａｎｄに基づいたしきい値は、高解像度化処理の周波数特性λｆ（ｕ，ｖ）に対して適切に設定される。したがって、該しきい値をしきい値処理に適用することによって、高解像度化処理によって増幅されたノイズを効果的に低減することができる。
【００５０】
さらに、しきい値処理によるノイズ低減処理を行って取得した画像を用いてノイズ低減フィルタを作成し、ウェーブレット空間上で適用する手法が知られている。
【００５１】
ノイズ低減フィルタλｓｕｂｂａｎｄ（ｘ，ｙ）は、しきい値処理によるノイズ低減処理を行って取得した画像をウェーブレット変換して得られるサブバンド画像ｗｓｕｂｂａｎｄ（ｘ，ｙ）と、上記してきたノイズ量σｓｕｂｂａｎｄに基づいて（１６）式により作成される。
【００５２】
【数１０】
2015166979.tif 000012
【００５３】
上記ノイズ低減フィルタはしきい値処理と同様に、高解像度画像をウェーブレット変換により分解したＬＬ以外のサブバンド画像に対して適用することで、高解像度化処理によって増幅されたノイズを低減することができる。
【００５４】
さらに、ノイズ低減処理を行ったサブバンド画像をウェーブレット逆変換により再合成することで、ノイズを低減した画像を生成する。このノイズ低減画像は、高解像度化処理の周波数特性に基づいて、高解像度化処理により増幅されたノイズが効果的に低減されたものである。
【００５５】
各実施例では、高解像度画像を以上説明したウェーブレット変換により周波数分解して得られた複数のサブバンド画像における各サブバンド画像にて分割された領域ごとにノイズ量を推定し、ノイズ低減処理を行う。そして、ノイズが低減された複数のサブバンド画像を再合成することで、高解像度化処理により増幅されたノイズが低減された画像を生成する。
【実施例１】
【００５６】
図１には、本発明の実施例１である画像処理方法の流れを示している。この画像処理方法は、画像処理装置に含まれるＣＰＵ等により構成されるコンピュータがコンピュータプログラムとして画像処理プログラムにしたがって実行する。
【００５７】
本実施例の画像処理方法は、撮像光学系を用いた撮像により生成された入力画像に対して高解像度化処理を行い、その後、高解像度化処理の周波数特性に基づいてノイズ低減処理を行う。本実施例では、高解像度化処理として、画像回復処理を行う場合について説明する。
【００５８】
ステップＳ１０１おいて、コンピュータは、カメラ等の撮像装置から撮影画像を入力画像として取得する。撮像装置からの入力画像の取得は、該撮像装置と画像処理装置とを有線または無線による通信を介して行ってもよいし、半導体メモリや光ディスク等の記憶媒体を介して行ってもよい。入力画像には、撮像時における撮像光学系の焦点距離、Ｆ値および撮像距離（フォーカスレンズ位置）と、撮像装置に設けられた撮像素子の画素サイズ等の撮像条件を示す情報（以下、撮像条件情報という）が付加されている。
【００５９】
図２は、入力画像を水平方向にＭ分割し、垂直方向にＮ分割したときの複数（Ｍ×Ｎ個）のシフトインバリアント領域を示している。撮像光学系の光学特性は画像（撮像面）全体において一様ではない。しかし、画像中の狭い領域であればその中の光学特性は一様であるとみなすことができる。本実施例では、この狭い領域をシフトインバリアント(shift-invariant)領域といい、シフトインバリアント領域（ｍ，ｎ）ごとに処理を行う。なお、シフトインバリアント領域は、予め撮像光学系の光学特性を考慮して決定しておくことが好ましい。分割数を増やしてシフトインバリアント領域を細かくすると、後述する画像回復処理等の精度は向上するが、一方で処理速度やデータ量は増えることになる。したがって、処理精度、処理速度およびデータ量のバランスを考慮してシフトバリアント領域の領域数を決定するとよい。
【００６０】
次に、ステップＳ１０２において、コンピュータは、入力画像に含まれるノイズ量（第１のノイズ量）σ０を算出する。入力画像に含まれるノイズ量σ０は、入力画像をウェーブレット変換して取得した複数のサブバンド画像のうち最も高周波側のサブバンド画像ＨＨ１におけるＭＡＤを（１４）式により求めることで算出する。また、サブバンド画像ＨＨ１におけるシフトインバリアント領域ごとにＭＡＤを求めることで、シフトインバリアント領域（ｍ，ｎ）ごとに入力画像のノイズ量を算出してもよい。また、撮像装置のＩＳＯ感度や輝度レベル等の撮像条件に対するノイズ量を予め測定して得られたノイズ量データをメモリに記憶させ、入力画像の撮像条件情報とシフトインバリアント領域の平均輝度値に合致するノイズ量を入力画像のノイズ量として取得してもよい。
【００６１】
なお、ステップＳ１０２は、ステップＳ１０５を実行する前であれば、ステップＳ１０３およびステップＳ１０４を実行した後に実行してもよい。
【００６２】
次に、ステップＳ１０３において、コンピュータは、撮像光学系の光学特性である光学伝達関数を取得する。光学伝達関数とは、Ｈ（ｕ，ｖ）で表されるように、２次元空間におけるＰＳＦの周波数特性である。すなわち、光学伝達関数はＰＳＦによってその特性が変化する。光学伝達関数は撮像光学系の状態によって膨大な量のバリエーションがある。例えば、焦点距離が可変であるズームレンズにおいては、焦点距離、Ｆ値、撮像距離、撮像素子（撮像面）上での位置、色（波長）等の組み合わせ分だけ光学伝達関数が存在する。さらに、撮像素子の画素サイズ等の撮像素子の特性と組み合わせた場合にもバリエーションがある。したがって、現実的には、多数の組み合わせの中から画像回復処理に必要なものだけを保持しておけばよい。なお、保持している撮像条件以外の撮像条件での光学伝達関数は、補間処理によって生成してもよい。これら光学伝達関数データは予め光学ＣＡＤ等で計算しておき、画像回復処理に適したデータベースにしてメモリに保持しておく。光学伝達関数データベースの撮像面上での位置は、図２に示した各シフトインバリアント領域に対応するようにしておくことが望ましい。また、光学伝達関数は測定によって求め、データとして保持しておき、取得してもよい。
【００６３】
そして、取得した入力画像から撮像時の焦点距離、Ｆ値、撮像距離、画素サイズ等の情報を取得し、予め保持された光学伝達関数データベースの中から上記に合致した光学伝達関数を取得する。
【００６４】
次に、ステップＳ１０４において、コンピュータは、高解像度化処理の周波数特性を算出する。ユーザは画像回復処理の強さを指定し、取得した光学伝達関数を用いて（４）式により高解像度化処理の周波数特性λｆ（ｕ，ｖ）を算出する。
【００６５】
図２に示すシフトインバリアント領域（ｍ，ｎ）を代表する光学伝達関数を、
Ｈｍ，ｎ（ｕ，ｖ）
とする。このとき、（４）式は、
【００６６】
【数１１】
2015166979.tif 000013
【００６７】
と書き換えることができる。
【００６８】
（１７）式により、シフトインバリアント領域（ｍ，ｎ）ごとに高解像度化処理の周波数特性λｆｍ，ｎ（ｕ，ｖ）を算出することができる。
【００６９】
画像回復処理の強さはユーザが分かりやすい０〜１００等のような数値範囲でよい。このとき、ユーザが指定した値に対して、実際に適用する画像回復処理の信号増幅量を設定する。この信号増幅量は、（１７）式の絶対値（｜λｆｍ，ｎ（ｕ，ｖ）｜）として計算され、Γｍ，ｎによって制御を行う。
【００７０】
信号増幅量は、例えば特定の周波数位置（Ｕ，Ｖ）における値｜λｍ，ｎｆ（Ｕ，Ｖ）｜としてもよいし、｜λｆｍ，ｎ（ｕ，ｖ）｜の帯域内の最大値や平均値としてもよい。
【００７１】
次に、ステップＳ１０５において、コンピュータは、入力画像に対して、算出した高解像度化処理の周波数特性λｆｍ，ｎ（ｕ，ｖ）に基づいて画像回復処理を行い、高解像度画像（第１の画像）を取得（生成）する。このとき、シフトインバリアント領域（ｍ，ｎ）ごとに算出した高解像度化処理の周波数特性λｆｍ，ｎ（ｕ，ｖ）を用いることによって、シフトインバリアント領域（ｍ，ｎ）ごとに最適な画像回復処理を行うことが可能となる。さらに言えば、シフトインバリアント領域の切り替わりでブロックノイズ等が発生しないような補間処理を用いると、より良好な画像回復処理による高解像度化処理を行うことが可能となる。また、画像回復処理は、（１７）式をフーリエ逆変換して得られる実空間のフィルタを用いて、入力画像に畳み込む（コンボリューションする）ことで行ってもよい。この場合、予め生成された実空間のフィルタを光学伝達関数データベースに格納して用いることも可能である。
【００７２】
次に、ステップＳ１０６において、コンピュータは、高解像度画像に対してウェーブレット変換（離散ウェーブレット変換）を行う。これにより、高解像度画像は、図３に示すような互いに異なる複数の周波数帯域成分ＨＨ，ＨＬ，ＬＨ，ＬＬとしての複数のサブバンド画像（第２の画像）に周波数分解される。ウェーブレット変換で用いる基底関数は予め数列データとしてメモリに記憶しておき、選択した基底関数の数列を取得してウェーブレット変換を行う。本実施例では、変換レベルを３としたときの処理の例について説明するが、変換レベルは３に限定されるものではなく、変更が可能である。
【００７３】
ウェーブレット変換により分解して得られた複数のサブバンド画像は、図４に示すように、画像回復処理時に入力画像を分割して処理したときと同様に、垂直方向にＭ分割され、水平方向にＮ分割される。サブバンド画像における各分割領域は、撮像光学系の光学特性（光学伝達関数：ＯＴＦ）に基づいて分割された領域に相当する。そして、各サブバンド画像おけるＭ×Ｎ個の分割領域における分割領域（ｍ，ｎ）ごとに処理を行う。
【００７４】
次に、ステップＳ１０７において、コンピュータは、入力画像のシフトインバリアント領域（ｍ，ｎ）で行った高解像度化処理の周波数特性λｆｍ，ｎ（ｕ，ｖ）を取得し、サブバンド画像の分割領域（ｍ，ｎ）における増幅量Ａｓｕｂｂａｎｄ（ｍ，ｎ）を算出する。光学伝達関数データベースに実空間のフィルタを格納した場合の増幅量は、予め実空間のフィルタとともに格納したものを用いてもよいし、実空間のフィルタをフーリエ変換して求めてもよい。
【００７５】
上述したように撮像光学系の光学伝達関数Ｈｍ，ｎ（ｕ，ｖ）はシフトインバリアント領域ごとに異なるため、高解像度化処理は入力画像におけるシフトインバリアント領域ごとの周波数特性λｆｍ，ｎ（ｕ，ｖ）に基づいて行われる。このため、高解像度画像には、シフトインバリアント領域における高解像度化処理の周波数特性λｆｍ，ｎ（ｕ，ｖ）に基づいて増幅されたノイズを含んでいる。つまり、高解像度化処理によるノイズの増幅量は、各サブバンド画像の複数の分割領域において一様ではなく、分割領域ごとに異なる。
【００７６】
したがって、サブバンド画像の分割領域（ｍ，ｎ）における増幅量Ａｓｕｂｂａｎｄ（ｍ，ｎ）を算出する。このとき、（１３）式は、
【００７７】
【数１２】
2015166979.tif 000014
【００７８】
と書き換えることができる。
【００７９】
図５には、サブバンド画像の分割領域（ｍ，ｎ）における増幅量Ａｓｕｂｂａｎｄ（ｍ，ｎ）の算出結果の例を示している。
【００８０】
なお、周波数分解の周波数特性Ψｓｕｂｂａｎｄ（ｕ，ｖ）は、ウェーブレット変換時に取得した基底関数の数列をフーリエ変換することで取得する。また、予めウェーブレット変換の各基底関数の数列からサブバンドの周波数特性Ψｓｕｂｂａｎｄ（ｕ，ｖ）を算出して、メモリに記憶しておき、これを取得してもよい。
【００８１】
また、後述するノイズ低減処理（ステップＳ１０９）はＬＬ以外のサブバンド画像に対して行われるため、ＬＬにおけるノイズの増幅量の算出（ステップＳ１０７）と後述するノイズ量の算出（ステップＳ１０８）は省略することができる。
【００８２】
次に、ステップＳ１０８において、コンピュータは、サブバンド画像の各分割領域（ｍ，ｎ）におけるノイズ量（第２のノイズ量）を取得する。ここにいうノイズ量とは、ノイズ分布の標準偏差である。上述したように入力画像のシフトインバリアント領域（ｍ，ｎ）ごとに行った高解像度化処理の周波数特性λｆｍ，ｎ（ｕ，ｖ）は、シフトインバリアント領域ごとの撮像光学系の光学伝達関数Ｈｍ，ｎ（ｕ，ｖ）に基づいて算出される。このため、サブバンド画像の分割領域ごとにノイズ量が異なる。よって、高解像度化処理の周波数特性に基づいて、サブバンド画像の分割領域（ｍ，ｎ）ごとのノイズ量σｓｕｂｂａｎｄ（ｍ，ｎ）を取得する必要がある。
【００８３】
サブバンド画像の分割領域（ｍ，ｎ）ごとのノイズ量σｓｕｂｂａｎｄ（ｍ，ｎ）は、ステップＳ１０２で取得した入力画像に含まれるノイズ量σ０と、ステップＳ１０７で算出したサブバンド画像の分割領域ごとの増幅量Ａｓｕｂｂａｎｄ（ｍ，ｎ）とを用いて算出する。このとき、（１２）式は、
【００８４】
【数１３】
2015166979.tif 000015
【００８５】
と書き換えることができる。
【００８６】
次に、ステップＳ１０９において、コンピュータは、サブバンド画像の分割領域（ｍ，ｎ）ごとに取得したノイズ量に基づいて、該分割領域ごとにノイズ低減処理を行う。本実施例では、ノイズ低減処理として、しきい値処理を用いる。
【００８７】
しきい値処理では、まずコンピュータは、しきい値を設定する。しきい値の設定は、重みパラメータとして（２ｌｏｇ（Ｑｓｕｂｂａｎｄ））１／２を用い、算出したノイズ量に基づいて設定する。Ｑｓｕｂｂａｎｄはサブバンド画像ごとのデータ数である。
【００８８】
次に、コンピュータは、ＬＬ以外のサブバンド画像において、取得したしきい値に基づいて、サブバンド画像の分割領域（ｍ，ｎ）内の位置（ｉ，ｊ）ごとのウェーブレット係数ｗｍ，ｎｓｕｂｂａｎｄ（ｉ，ｊ）のうち、しきい値以下の係数を０に置き換える。これにより、しきい値処理を行う。このとき、（１１）式は、
【００８９】
【数１４】
2015166979.tif 000016
【００９０】
と書き換えることができる。
【００９１】
設定したしきい値は、ステップＳ１０８で算出したサブバンド画像の分割領域（ｍ，ｎ）におけるノイズ量に基づくものである。つまり、しきい値処理によるノイズ低減処理は、入力画像のシフトインバリアント領域ごとの高解像度化処理によるノイズの増幅量に対して適切に設定されている。したがって、取得したしきい値を用いてしきい値処理を行うことで、シフトインバリアント領域ごとの高解像度化処理の周波数特性に基づいて、増幅ノイズを効果的に低減することができる。
【００９２】
次に、ステップＳ１１０において、コンピュータは、ノイズ低減処理後の複数のサブバンド画像ＨＨ，ＨＬ，ＬＨ，ＬＬをウェーブレット逆変換により再合成する。これにより、ステップＳ１１１において、ノイズ低減処理後の画像（第３の画像）を取得する。ノイズ低減処理後の画像は、シフトインバリアント領域ごとの高解像度化処理の周波数特性に基づいて、高解像度化処理により増幅されたノイズが効果的に低減されたものである。
【００９３】
なお、本実施例では、高解像度化処理として画像回復処理を行う場合について説明したが、高解像度化処理の周波数特性を取得できれば、エッジ強調処理や超解像処理等の高解像度化処理においても、本実施例と同様の画像処理方法を適用することができる。
【実施例２】
【００９４】
次に、本発明の実施例２について説明する。実施例１では、図１のステップＳ１０９におけるノイズ低減処理としてしきい値処理を用いたが、本実施例では、高解像度画像に対して、ノイズ低減フィルタを用いたノイズ低減処理を行う。本実施例における画像処理方法は、コンピュータが、図６のフローチャートに示した順序で画像処理プログラムに従って実行する。
【００９５】
ノイズ低減フィルタは、予めしきい値処理によるノイズ低減処理を行った画像を用いて作成する。したがって、まずコンピュータは、ステップＳ２０１〜Ｓ２１１において、シフトインバリアント領域ごとに高解像度化処理を行って取得した高解像度画像に対して、しきい値処理によるノイズ低減処理を行い、ノイズ低減処理後の画像を取得する。ここにいうしきい値処理とは、実施例１で詳細に説明したように、シフトインバリアント領域ごとの高解像度化処理の周波数特性に基づいて、増幅ノイズを効果的に低減したものである。
【００９６】
次に、ステップＳ２１２において、コンピュータは、しきい値処理によるノイズ低減処理を行って取得した画像に対してウェーブレット変換を行い、複数の周波数帯域成分ＨＨ，ＨＬ，ＬＨ，ＬＬとしての複数のサブバンド画像に周波数分解する。このとき、ウェーブレット変換の基底関数は後述するステップＳ２１３で用いるものを用いる。なお、ステップＳ２１２は、ステップＳ２１１を実行した後からステップＳ２１６を実行する前であればいつ行ってもよい。また、ステップＳ２０６とステップＳ２１３で用いる基底関数が同じである場合にはステップＳ２１０〜Ｓ２１２およびステップＳ２１４、Ｓ２１５を省略してもよい。
次に、ステップＳ２１３において、コンピュータは、高解像度画像に対してウェーブレット変換を行い、複数の周波数帯域成分ＨＨ，ＨＬ，ＬＨ，ＬＬとしての複数のサブバンド画像に周波数分解する。
次に、ステップＳ２１４〜Ｓ２１５において、コンピュータは、高解像度化処理の周波数特性に基づいて、ステップＳ２１３で取得したサブバンド画像の分割領域（ｍ，ｎ）における増幅量Ａｓｕｂｂａｎｄ（ｍ，ｎ）を算出し、ノイズ量σｓｕｂｂａｎｄ（ｍ，ｎ）を取得する。本実施例における図６のステップＳ２１４〜Ｓ２１５は実施例１における図１のステップＳ１０７〜Ｓ１０８と同様の処理であるため、詳細な説明は省略する。
【００９７】
次に、ステップＳ２１６では、コンピュータは、ノイズ低減フィルタλｍ,ｎｓｕｂｂａｎｄ（ｉ，ｊ）を作成する。具体的には、サブバンド画像のウェーブレット係数ｗｍ,ｎｓｕｂｂａｎｄ（ｉ，ｊ）と、ステップＳ２１５で算出したサブバンド画像の分割領域（ｍ，ｎ）におけるノイズ量σｓｕｂｂａｎｄ（ｍ，ｎ）とを用いて作成する。このとき、（１６）式は、
【００９８】
【数１５】
2015166979.tif 000017
【００９９】
と書き代えることができる。この際、ステップＳ２０６とステップＳ２１３で用いる基底関数が同じである場合には、ステップＳ２０８で算出したノイズ量を用いてノイズ低減フィルタを作成してもよい。なお、ノイズ低減処理はＬＬ以外のサブバンド画像に対して行うため、ＬＬにおけるノイズ低減フィルタの作成は省略することができる。
【０１００】
次に、ステップＳ２１７において、コンピュータは、上記の手法により作成したノイズ低減フィルタを、ステップＳ２１３で取得したＬＬ以外のサブバンド画像の分割領域ごとに適用し、ノイズ低減処理を行う。
次に、ステップＳ２１８において、コンピュータは、ノイズ低減処理後の複数のサブバンド画像ＨＨ，ＨＬ，ＬＨ，ＬＬをウェーブレット逆変換により再合成する。これにより、ステップＳ２１９において、コンピュータは、ノイズ低減処理後の画像を取得する。
【０１０１】
作成したノイズ低減フィルタは、算出したサブバンド画像の分割領域（ｍ，ｎ）におけるノイズ量に基づいたものである。つまり、上記したノイズ低減フィルタによるノイズ低減処理は、入力画像のシフトインバリアント領域ごとの高解像度化処理によるノイズの増幅量に対して適切に設定されている。したがって、作成したノイズ低減フィルタを用いてノイズ低減処理を行うことで、シフトインバリアント領域ごとの高解像度化処理の周波数特性に基づいて、高解像度化処理により増幅されたノイズを効果的に低減することができる。
【実施例３】
【０１０２】
次に、本発明の実施例３として、上述した実施例１または実施例２の画像処理方法（画像処理プログラム）を実行する画像処理装置（コンピュータ）を内蔵した撮像装置について説明する。図７には、該撮像装置の構成を示している。
【０１０３】
撮像光学系２０１により形成された被写体の光学像を撮像する撮像素子２０２は、該光学像を電気信号にして変換してアナログ信号を出力する。撮像素子２０２から出力されたアナログ信号は、Ａ／Ｄ変換器２０３によりデジタル信号に変換される。該デジタル信号は、画像処理部２０４に入力される。
【０１０４】
画像処理部２０４は、入力されたデジタル信号に対して各種画像処理を行うことにより撮影画像（入力画像）を生成する。
【０１０５】
また、画像処理部２０４は、表示部２０５に、ユーザが高解像度化処理およびノイズ低減処理を行うか否かを決定するためのＧＵＩと、高解像度化処理の強さを決定するためのＧＵＩとを表示する。ＧＵＩにおいてユーザが高解像度化処理およびノイズ低減処理を実行すると決定し、さらに高解像度化処理の強さを決定すると、画像処理部２０４は、入力画像に対して実施例１または実施例２で説明した高解像度化処理およびノイズ低減処理を実行する。また、予めユーザが高解像度化処理およびノイズ低減処理を実行すると決定するとともに高解像度化処理の強さを決定しておき、取得した入力画像に対して自動的に高解像度化処理およびノイズ低減処理を実行するようにしてもよい。
【０１０６】
画像処理部２０４に高解像度化処理を行わせるために、状態検知部２０７は、入力画像を生成する際における撮像条件情報を得る。状態検知部２０７は、撮像条件情報をシステムコントローラ２１０から得てもよい。また、撮像光学系に関する撮像条件情報については、撮像光学系２０１の絞り２０１ａの動作やズームレンズ（図示せず）およびフォーカスレンズ２０１ｂの移動を制御する撮像光学系制御部２０６から得てもよい。撮像条件情報は入力画像に付加され、該入力画像は画像処理部２０４に接続されている画像記録媒体２０９に記録される。
【０１０７】
記憶部２０８は、撮像光学系２０１の焦点距離、Ｆ値および撮像距離と、撮像素子２０２（撮像面）上での位置と、撮像素子２０２の画素サイズ等の撮像条件の組み合わせごとに光学伝達関数の情報を保持している。
【０１０８】
画像処理部２０４は、高解像度化処理およびノイズ低減処理により生成した出力画像（第３の画像）を半導体メモリや光ディスク等の画像記録媒体２０９に出力して記録させたり、表示部２０５に出力して表示させたりする。
【０１０９】
以上説明した一連の動作は、ＣＰＵやＭＰＵ等によって構成されるシステムコントローラ２１０により制御される。
【０１１０】
以上説明した各実施例は代表的な例にすぎず、本発明の実施に際しては、各実施例に対して種々の変形や変更が可能である。
【産業上の利用可能性】
【０１１１】
ノイズが少ない高解像度の画像を生成可能な画像処理装置や撮像装置を提供できる。
【符号の説明】
【０１１２】
２０１  撮像光学系
２０２  撮像素子
２０４  画像処理部
【図１】
2015166979.tif 000018
【図２】
2015166979.tif 000019
【図３】
2015166979.tif 000020
【図４】
2015166979.tif 000021
【図５】
2015166979.tif 000022
【図６】
2015166979.tif 000023
【図７】
2015166979.tif 000024