(57)【要約】      （修正有）【課題】画像処理用のTotalVariationフィルタにおいて、繰り返し演算に要する時間を短縮する。【解決手段】入力画像から関数Ｆ（ｕ）＝J（ｕ）＋λΣ（ｕi,j−ｆi,j）２を最小とするような骨格成分を求めるTotalVariationフィルタにおいて、Ｊ（ｕ）としてＪ（ｕ）＝Σ（｜ｕi+1,j−ｕi,j｜＋｜ｕi,j+1−ｕi,j｜＋｜ｕi+1,j+1−ｕi,j｜＋｜ｕi,j+1−ｕi+1,j｜）を用いる。この場合、関数Ｆ（ｕ）＝J（ｕ）＋λΣ（ｕi,j−ｆi,j）２を最小化するアルゴリズムにおいて、注目画素値ｕi,jに隣接する８画素値とｕi,jとの差分の符号ｓｇｎ（）を取ったものを８個加算して修正項Ａを求める。これにより、繰り返し演算の一回の計算に用いる情報量が２倍になるため、収束速度が２倍に向上する。【選択図】図６
2011180717.tif 000002
【特許請求の範囲】
【請求項１】
入力画像から関数Ｆ（ｕ）＝J（ｕ）＋λΣ（ｕi,j−ｆi,j）２を最小とするような骨格成分を求める画像処理用のTotal Variation フィルタにおいて、
前記Ｊ（ｕ）としてＪ（ｕ）＝Σ（｜ｕi+1,j−ｕi,j｜＋｜ｕi,j+1−ｕi,j｜＋｜ｕi+1,j+1−ｕi,j｜＋｜ｕi,j+1−ｕi+1,j｜）を用いることを特徴とする画像処理用のTotal Variation フィルタ。
但し、ｆi,jは入力画素値、ｕi,jは注目画素値、ｕi+1,j、ｕi,j+1、ｕi+1,j+1は隣接画素値である。
【請求項２】
入力画像から関数Ｆ（ｕ）＝J（ｕ）＋λΣ（ｕi,j−ｆi,j）２を最小とするような骨格成分を求める画像処理用のTotal Variation フィルタにおいて、
前記Ｊ（ｕ）としてＪ（ｕ）＝Σ（｜ｕi+1,j−ｕi,j｜＋｜ｕi,j+1−ｕi,j｜＋｜ｕi+1,j+1−ｕi,j｜）を用いることを特徴とする画像処理用のTotal Variation フィルタ。
但し、ｆi,jは入力画素値、ｕi,jは注目画素値、ｕi+1,j、ｕi,j+1、ｕi+1,j+1は隣接画素値である。
【請求項３】
入力画像から関数Ｆ（ｕ）＝J（ｕ）＋λΣ（ｕi,j−ｆi,j）２を最小とするような骨格成分を求める画像処理用のTotal Variation フィルタにおいて、
前記Ｊ（ｕ）としてＪ（ｕ）＝Σ（｜ｕi+1,j−ｕi,j｜＋｜ｕi,j+1−ｕi,j｜＋｜ｕi,j+1−ｕi+1,j｜）を用いることを特徴とする画像処理用のTotal Variation フィルタ。
但し、ｆi,jは入力画素値、ｕi,jは注目画素値、ｕi+1,j、ｕi,j+1は隣接画素値である。
【請求項４】
入力画像から画像の骨格成分を求める画像処理用のTotal Variation フィルタにおいて、
Ａ＝sgn（ｕi,j−ｕi+1,j）＋sgn（ｕi,j−ｕi,j+1）＋sgn（ｕi,j−ｕi+1,j+1）＋sgn（ｕi,j−ｕi-1,j）＋sgn（ｕi,j−ｕi-1,j-1）＋sgn（ｕi,j−ｕi,j-1）＋sgn（ｕi,j−ｕi-1,j+1）＋sgn（ｕi,j−ｕi+1,j-1）を求める第１の手段と、
現在のｕi,jに−α｛２λ（ｕi,j−ｆi,j）＋Ａ｝を加算して、ｕi,jを更新する第２の手段と、
前記第１の手段と前記第２の手段の処理をＮ回繰り返させることにより、骨格成分を得る第３の手段と、を備えたことを特徴とする画像処理用のTotal Variation フィルタ。
但し、ｆi,jは入力画素値、ｕi,jは注目画素値、ｕi+1,j、ｕi,j+1、ｕi+1,j+1は隣接画素値、λ、αはそれぞれ正の定数である。
【請求項５】
入力画像から画像の骨格成分を求める画像処理用のTotal Variation フィルタにおいて、
Ａ＝sgn（ｕi,j−ｕi+1,j）＋sgn（ｕi,j−ｕi,j+1）＋sgn（ｕi,j−ｕi+1,j+1）＋sgn（ｕi,j−ｕi-1,j）＋sgn（ｕi,j−ｕi-1,j-1）＋sgn（ｕi,j−ｕi,j-1）を求める第１の手段と、
現在のｕi,jに−α｛２λ（ｕi,j−ｆi,j）＋Ａ｝を加算して、ｕi,jを更新する第２の手段と、
前記第１の手段と前記第２の手段の処理をＮ回繰り返させることにより、骨格成分を得る第３の手段と、を備えたことを特徴とする画像処理用のTotal Variation フィルタ。
但し、ｆi,jは入力画素値、ｕi,jは注目画素値、ｕi+1,j、ｕi,j+1、ｕi+1,j+1は隣接画素値、λ、αはそれぞれ正の定数である。
【請求項６】
入力画像から画像の骨格成分を求める画像処理用のTotal Variation フィルタにおいて、
Ａ＝sgn（ｕi,j−ｕi+1,j）＋sgn（ｕi,j−ｕi,j+1）＋sgn（ｕi,j−ｕi-1,j）＋sgn（ｕi,j−ｕi,j-1）＋sgn（ｕi,j−ｕi-1,j+1）＋sgn（ｕi,j−ｕi+1,j-1）を求める第１の手段と、
現在のｕi,jに−α｛２λ（ｕi,j−ｆi,j）＋Ａ｝を加算して、ｕi,jを更新する第２の手段と、
前記第１の手段と前記第２の手段の処理をＮ回繰り返させることにより、骨格成分を得る第３の手段と、を備えたことを特徴とする画像処理用のTotal Variation フィルタ。
但し、ｆi,jは入力画素値、ｕi,jは注目画素値、ｕi+1,j、ｕi,j+1は隣接画素値、λ、αはそれぞれ正の定数である。
【請求項７】
入力画像から画像の骨格成分を求めるための画像処理プログラムであって、
コンピュータを、
Ａ＝sgn（ｕi,j−ｕi+1,j）＋sgn（ｕi,j−ｕi,j+1）＋sgn（ｕi,j−ｕi+1,j+1）＋sgn（ｕi,j−ｕi-1,j）＋sgn（ｕi,j−ｕi-1,j-1）＋sgn（ｕi,j−ｕi,j-1）＋sgn（ｕi,j−ｕi-1,j+1）＋sgn（ｕi,j−ｕi+1,j-1）を求める第１の手段と、
現在のｕi,jに−α｛２λ（ｕi,j−ｆi,j）＋Ａ｝を加算して、ｕi,jを更新する第２の手段と、
前記第１の手段と前記第２の手段の処理をＮ回繰り返させることにより、骨格成分を得る第３の手段として機能させることを特徴とする画像処理プログラム。
但し、ｆi,jは入力画素値、ｕi,jは注目画素値、ｕi+1,j、ｕi,j+1、ｕi+1,j+1は隣接画素値、λ、αはそれぞれ正の定数である。
【請求項８】
入力画像から画像の骨格成分を求めるための画像処理プログラムであって、
コンピュータを、
Ａ＝sgn（ｕi,j−ｕi+1,j）＋sgn（ｕi,j−ｕi,j+1）＋sgn（ｕi,j−ｕi+1,j+1）＋sgn（ｕi,j−ｕi-1,j）＋sgn（ｕi,j−ｕi-1,j-1）＋sgn（ｕi,j−ｕi,j-1）を求める第１の手段と、
現在のｕi,jに−α｛２λ（ｕi,j−ｆi,j）＋Ａ｝を加算して、ｕi,jを更新する第２の手段と、
前記第１の手段と前記第２の手段の処理をＮ回繰り返させることにより、骨格成分を得る第３の手段として機能させることを特徴とする画像処理プログラム。
但し、ｆi,jは入力画素値、ｕi,jは注目画素値、ｕi+1,j、ｕi,j+1、ｕi+1,j+1は隣接画素値、λ、αはそれぞれ正の定数である。
【請求項９】
入力画像から画像の骨格成分を求めるための画像処理プログラムであって、
コンピュータを、
Ａ＝sgn（ｕi,j−ｕi+1,j）＋sgn（ｕi,j−ｕi,j+1）＋sgn（ｕi,j−ｕi-1,j）＋sgn（ｕi,j−ｕi,j-1）＋sgn（ｕi,j−ｕi-1,j+1）＋sgn（ｕi,j−ｕi+1,j-1）を求める第１の手段と、
現在のｕi,jに−α｛２λ（ｕi,j−ｆi,j）＋Ａ｝を加算して、ｕi,jを更新する第２の手段と、
前記第１の手段と前記第２の手段の処理をＮ回繰り返させることにより、骨格成分を得る第３の手段として機能させることを特徴とする画像処理プログラム。
但し、ｆi,jは入力画素値、ｕi,jは注目画素値、ｕi+1,j、ｕi,j+1は隣接画素値、λ、αはそれぞれ正の定数である。
【発明の詳細な説明】
【技術分野】
【０００１】
本発明は、画像処理用のTotal Variation フィルタおよび画像処理プログラムに関し、特にその高速化を行って繰り返し演算に要する時間を短縮するものに関する。
【背景技術】
【０００２】
画像信号処理において、Total Variation フィルタは、ノイズの除去や超解像画像拡大などで非常に有効な手段となっている。Total Variation フィルタは、通常のフィルタと異なり、非線形の繰り返し演算によって、画像信号処理を行う（特許文献１〜３、非特許文献１、２参照）。
【先行技術文献】
【特許文献】
【０００３】
【特許文献１】特開２００８−１７２４３１号公報
【特許文献２】特開２００９−５２５２号公報
【特許文献３】特開２００９−２０６０５号公報
【非特許文献】
【０００４】
【非特許文献１】L. Rudin, S. Osher, and E. Fetami, "Nonlinear total Variation Based Noise Removal Algorithm", Physica D, Vol.60, pp.259-268, 1992.
【非特許文献２】A. Beck and M. Teboulle, "Fast Gradient-Based Algorithms for Constrained Total Variation Image Denoising and Deblurring Problems", IEEE Trans. on Image Processing, Vol.18, No.11, pp.2419-2434, Nov.2009.
【発明の概要】
【発明が解決しようとする課題】
【０００５】
Total Variation フィルタは繰り返し演算（具体的には１００〜２００回）が必要であり、処理に非常に時間がかかる。このため、その優れた特性にもかかわらず静止画の処理のみに用いられており、テレビや映画などの動画の処理には用いることができなかった。これを動画処理に用いるためには、繰り返し演算の回数の削減（収束速度の高速化）が必須となる。
【０００６】
本発明は上記点に鑑みて、Total Variation フィルタにおいて、その高速化を行って繰り返し演算に要する時間を短縮することを目的とする。
【課題を解決するための手段】
【０００７】
上記目的を達成するため、発明者らは、注目画素の上下左右の隣接画素の差分の絶対値和に斜め方向の隣接画素の差分の絶対値和を付け加えれば、繰り返し演算の一回の計算に用いる情報量が２倍になるため、収束速度が２倍に向上する、ということを考えた。
【０００８】
そこで、請求項１に記載の発明は、
入力画像から関数Ｆ（ｕ）＝J（ｕ）＋λΣ（ｕi,j−ｆi,j）２を最小とするような骨格成分を求める画像処理用のTotal Variation フィルタにおいて、
前記Ｊ（ｕ）としてＪ（ｕ）＝Σ（｜ｕi+1,j−ｕi,j｜＋｜ｕi,j+1−ｕi,j｜＋｜ｕi+1,j+1−ｕi,j｜＋｜ｕi,j+1−ｕi+1,j｜）を用いることを特徴とする。
【０００９】
請求項４に記載の発明は、
入力画像から画像の骨格成分を求める画像処理用のTotal Variation フィルタにおいて、
Ａ＝sgn（ｕi,j−ｕi+1,j）＋sgn（ｕi,j−ｕi,j+1）＋sgn（ｕi,j−ｕi+1,j+1）＋sgn（ｕi,j−ｕi-1,j）＋sgn（ｕi,j−ｕi-1,j-1）＋sgn（ｕi,j−ｕi,j-1）＋sgn（ｕi,j−ｕi-1,j+1）＋sgn（ｕi,j−ｕi+1,j-1）を求める第１の手段と、
現在のｕi,jに−α｛２λ（ｕi,j−ｆi,j）＋Ａ｝を加算して、ｕi,jを更新する第２の手段と、
前記第１の手段と前記第２の手段の処理をＮ回繰り返させることにより、骨格成分を得る第３の手段と、を備えたことを特徴とする。
【００１０】
請求項７に記載の発明は、
入力画像から画像の骨格成分を求めるための画像処理プログラムであって、
コンピュータを、
Ａ＝sgn（ｕi,j−ｕi+1,j）＋sgn（ｕi,j−ｕi,j+1）＋sgn（ｕi,j−ｕi+1,j+1）＋sgn（ｕi,j−ｕi-1,j）＋sgn（ｕi,j−ｕi-1,j-1）＋sgn（ｕi,j−ｕi,j-1）＋sgn（ｕi,j−ｕi-1,j+1）＋sgn（ｕi,j−ｕi+1,j-1）を求める第１の手段と、
現在のｕi,jに−α｛２λ（ｕi,j−ｆi,j）＋Ａ｝を加算して、ｕi,jを更新する第２の手段と、
前記第１の手段と前記第２の手段の処理をＮ回繰り返させることにより、骨格成分を得る第３の手段として機能させることを特徴とするプログラム。
【００１１】
なお、上記において、ｆi,jは入力画素値、ｕi,jは注目画素値、ｕi+1,j、ｕi,j+1、ｕi+1,j+1は隣接画素値、λ、αはそれぞれ正の定数である。
【００１２】
また、請求項１に記載のＪ（ｕ）の代わりに、請求項２に記載の発明のように、Ｊ（ｕ）＝Σ（｜ｕi+1,j−ｕi,j｜＋｜ｕi,j+1−ｕi,j｜＋｜ｕi+1,j+1−ｕi,j｜）を用いてもよく、請求項３に記載の発明のように、Ｊ（ｕ）＝Σ（｜ｕi+1,j−ｕi,j｜＋｜ｕi,j+1−ｕi,j｜＋｜ｕi,j+1−ｕi+1,j｜）を用いてもよい。その場合、請求項４の構成の代わりに請求項５、６の構成とし、請求項７の構成の代わりに請求項８、９の構成とすることができる。
【図面の簡単な説明】
【００１３】
【図１】Total Variation フィルタの構成を示す図である。
【図２】従来法におけるＪ（ｕ）項を表す図である。
【図３】第１実施形態におけるＪ（ｕ）項を表す図である。
【図４】従来法における最急降下法において修正項Ａの求め方を説明するための説明図である。
【図５】本実施形態における最急降下法において修正項Ａの求め方を説明するための説明図である。
【図６】画像処理プログラムを示すフローチャートである。
【図７】第１実施形態の効果を示す図である。
【図８】第２実施形態におけるＪ（ｕ）項を表す図である。
【図９】第２実施形態における最急降下法において修正項Ａの求め方を説明するための説明図である。
【図１０】第３実施形態におけるＪ（ｕ）項を表す図である。
【図１１】第３実施形態における最急降下法において修正項Ａの求め方を説明するための説明図である。
【発明を実施するための形態】
【００１４】
（第１実施形態）
図１に、Total Variation フィルタの構成を示す。Total Variation フィルタは、入力画像から関数Ｆ（ｕ）＝J（ｕ）＋λΣ（ｕi,j−ｆi,j）２を最小とするような骨格成分を求めるために反復計算を行う。関数F(u)が最小となったときのｕは画像の骨格成分となる。この骨格成分を入力画像から単純に引くことによって微小振動であるテクスチャ成分が得られる。なお、λは適当な正の定数であり、その値を高く設定した場合、繰り返し演算後に得られる骨格成分は、原画像に近いものとなる。
【００１５】
評価関数J（ｕ）として、従来法では、
Ｊ（ｕ）＝Σ（｜ｕi+1,j−ｕi,j｜＋｜ｕi,j+1−ｕi,j｜    ・・・（１）
を用いていた。これに対し、本実施形態では、
Ｊ（ｕ）＝Σ（｜ｕi+1,j−ｕi,j｜＋｜ｕi,j+1−ｕi,j｜＋｜ｕi+1,j+1−ｕi,j｜＋｜ｕi,j+1−ｕi+1,j｜）    ・・・（２）
を用いる。
【００１６】
図２に従来法におけるＪ（ｕ）項を表す。注目画素値ｕｉ,ｊとその上下の画素値の差分をとり、その絶対値を式（１）で表わされるようにｉ，ｊを１から画素数Ｍまで加算して、Ｊ（ｕ）の値を得る。
【００１７】
図３に本実施形態におけるＪ（ｕ）項を表す。注目画素値ｕｉ,ｊとその上下の画素値の差分および斜め方向の２個の差分を取り、それらの絶対値を式（２）で表わされるようにｉ，ｊを１から画素数Ｍまで加算して、Ｊ（ｕ）の値を得る。
【００１８】
関数Ｆ（ｕ）＝J（ｕ）＋λΣ（ｕi,j−ｆi,j）２を最小化するアルゴリズムとしては、次のような最急降下法を用いる。
【００１９】
ｕi,j＝ｕi,j−α{２λ（ｕi,j−ｆi,j）＋Ａ}    ・・・（３）
Ａ＝sgn（ｕi,j−ｕi+1,j）＋sgn（ｕi,j−ｕi,j+1）＋sgn（ｕi,j−ｕi+1,j+1）＋sgn（ｕi,j−ｕi-1,j）＋sgn（ｕi,j−ｕi-1,j-1）＋sgn（ｕi,j−ｕi,j-1）＋sgn（ｕi,j−ｕi-1,j+1）＋sgn（ｕi,j−ｕi+1,j-1）      ・・・（４）
ただし、αは、ステップ幅で、適当な正の定数である。なお、式（３）において、左辺のｕｉ,ｊは更新後のものであり、右辺のｕｉ,ｊは現在のものである。
【００２０】
従来法では、上記した最急降下法において、図４に示すようにして修正項Ａを求める。この場合、注目画素値ｕｉ,ｊと上下左右の画素値との差分の符号ｓｇｎ（ ）を４個求める。
【００２１】
これに対し、本実施形態では、上記した最急降下法において、図５に示すようにして修正項Ａを求める。この場合、注目画素値ｕｉ,ｊの上下左右と斜め方向の８個の差分の符号ｓｇｎ（ ）を求め、式（４）のように、修正項Ａを求める。
【００２２】
Total Variation フィルタは、コンピュータを用いたソフトウェアにより実現することができる。図６にTotal Variation フィルタを用いた画像処理プログラムをフローチャートで示す。
【００２３】
ステップ１０１で演算回数Ｎが０に初期設定された後、ステップ１０２で修正項Ａが式（４）のように計算される。すなわち、注目画素値ｕi,jに隣接する８画素値とｕi,jとの差分の符号ｓｇｎ（ ）を取ったものを８個加算して値Ａを得る。ステップ１０３で画素値ｕi,jが−α{２λ（ｕi,j−ｆi,j）＋Ａ}によって新しい画素値ｕi,jに更新される。そして、ステップ１０４で演算回数Ｎがインクリメントされ、ステップ１０５で予め定められた値ＮstopにＮが達したか否かが判定される。Ｎが値Ｎstopに達していない場合はステップ１０２に戻る。Ｎが値Ｎstopに達した場合は画素値ｕi,jが最終骨格成分として出力され、またステップ１０６で入力画像ｆi,jからｕi,jが引き算されて、テクスチャ成分ｖi,jが出力される。
【００２４】
図７に本実施形態の効果を示す。縦軸は式（３）の最小化関数Ｆ（ｕ）の値、横軸は繰り返し回数Ｎの値を示す。本実施形態（図中、ＴＶＤ methodと記す）は従来法（図中、ＴＶＡ  methodと記す）にくらべて、半分の回数で収束している、すなわち収束速度が２倍になっていることが分かる。
【００２５】
以上述べたように、本実施形態によれば、画像処理において非常に高い性能を有する、Total Variation フィルタの演算の高速化を実現することができる。
【００２６】
（第２実施形態）
上記した実施形態では、評価関数J（ｕ）として、Ｊ（ｕ）＝Σ（｜ｕi+1,j−ｕi,j｜＋｜ｕi,j+1−ｕi,j｜＋｜ｕi+1,j+1−ｕi,j｜＋｜ｕi,j+1−ｕi+1,j｜）を用いるものを示したが、Ｊ（ｕ）＝Σ（｜ｕi+1,j−ｕi,j｜＋｜ｕi,j+1−ｕi,j｜＋｜ｕi+1,j+1−ｕi,j｜）を用いるようにしてもよい。この場合、図８に示すように、注目画素値ｕｉ,ｊとその上下の画素値の差分および斜め方向の１個の差分を取り、それらの絶対値を上記の式で表わされるようにｉ，ｊを１から画素数Ｍまで加算して、Ｊ（ｕ）の値を得る。また、この場合、最急降下法において、図９に示すように、注目画素値ｕｉ,ｊの上下左右と斜め方向の６個の差分の符号ｓｇｎ（ ）により修正項Ａを求める。具体的には、次のようにして修正項Ａを求める。
【００２７】
Ａ＝sgn（ｕi,j−ｕi+1,j）＋sgn（ｕi,j−ｕi,j+1）＋sgn（ｕi,j−ｕi+1,j+1）＋sgn（ｕi,j−ｕi-1,j）＋sgn（ｕi,j−ｕi-1,j-1）＋sgn（ｕi,j−ｕi,j-1）
（第３実施形態）
また、評価関数J（ｕ）として、Ｊ（ｕ）＝Σ（｜ｕi+1,j−ｕi,j｜＋｜ｕi,j+1−ｕi,j｜＋｜ｕi,j+1−ｕi+1,j｜）を用いるようにしてもよい。この場合、図１０に示すように、注目画素値ｕｉ,ｊとその上下の画素値の差分および斜め方向の１個の差分を取り、それらの絶対値を上記の式で表わされるようにｉ，ｊを１から画素数Ｍまで加算して、Ｊ（ｕ）の値を得る。また、この場合、最急降下法において、図１１に示すように、注目画素値ｕｉ,ｊの上下左右と斜め方向の６個の差分の符号ｓｇｎ（ ）により修正項Ａを求める。具体的には、次のようにして修正項Ａを求める。
【００２８】
Ａ＝sgn（ｕi,j−ｕi+1,j）＋sgn（ｕi,j−ｕi,j+1）＋sgn（ｕi,j−ｕi-1,j）＋sgn（ｕi,j−ｕi,j-1）＋sgn（ｕi,j−ｕi-1,j+1）＋sgn（ｕi,j−ｕi+1,j-1）
上記した第１実施形態では、画素値のＸ方向とＹ方向の隣接画素の差分の絶対値和に斜め方向の隣接画素の差分の絶対値和が２つ加わるの対し、第２、第３実施形態では、画素値のＸ方向とＹ方向の隣接画素の差分の絶対値和に斜め方向の隣接画素の差分の絶対値和が１つ加わるが、この第２、第３実施形態においても、従来法に比べれば、繰り返し演算の一回の計算に用いる情報量が多くなるため、収束速度を向上させることができる。
【図１】
2011180717.tif 000003
【図２】
2011180717.tif 000004
【図３】
2011180717.tif 000005
【図４】
2011180717.tif 000006
【図５】
2011180717.tif 000007
【図６】
2011180717.tif 000008
【図７】
2011180717.tif 000009
【図８】
2011180717.tif 000010
【図９】
2011180717.tif 000011
【図１０】
2011180717.tif 000012
【図１１】
2011180717.tif 000013