JPH04101282A - 空間周波数フィルタリング装置 - Google Patents

Abstract

(57)【要約】本公報は電子出願前の出願データであるた
め要約のデータは記録されません。

Description

【発明の詳細な説明】 ［産業上の利用分野コ 本発明は画像情報について強調、平滑化、劣化復元など
の処理を行なうための空間周波数フィルタリング装置に
関する。

［従来の技術］ 画像情報について強調、平滑化、劣化復元なと柾々の処
理を行なう場合、空間周波数のフィルタリングは有用で
ある。画像のフィルタリングの方法には２種類ある。１
つは、フーリエ変換により画像を空間周波数領域へ変換
し、各周波数ごとに適当な係数を掛けた後に逆フーリエ
変換して、画像に戻す方法である。以下、この方法をフ
ーリエ変換法と呼ぶことにする。他の１つは画像の注目
点（画Ｘ）の局所近傍領域について設定した係数マスク
と積和演算を行ない、値を置き換えるという操作を画像
の全面にわたって行なう方法で、いわゆる畳み込み演算
（コンボリユーシヨン）を行なう方法である。以下、こ
の方法をコンボリューション法と呼ぶことにする。なお
フーリエ変換の重要な性質としてコンボリューション定
理が知られており、上記両方法は空間周波数に対するフ
ィルタリングとして等価である。

フーリエ変換法を前提とした場合、空間周波数のフィル
ターは、基本的に画像と同じサイズで設計できる。した
がって自由度が大きく、様々な応用が考えられる。例え
ば下記文献１）では、内視鏡画像に対して輝度成分のフ
ィルタリングを行なうことにより、照明ムラを減らし、
観察中の体内組織の構造を強調し、同時に高周波ノイズ
の軽減を行なう方法が提案されている。また下記の文献
２）には劣化画像の復元法としてウィーナフィルターや
擬似逆行列による応用が紹介されているが、いずれもフ
ーリエ変換法を用いた処理である。

文献１　）　Ｎ、０ｈｙａａ＋、に、５ｕｚｕｋ１．Ｔ
、１ｌｏｎｄａ、Ｊ、Ｔｓｕｊｌｕｅｈｉ、Ｒ，Ｏｎｏ
　　ａｎｄ　　Ｓ、Ｉｋｃｄａ　　Ｏｐｔ　　Ｃｏｍ１
．．５５，２４２．（１９８５）。

文献２　）　Ｖ、に、Ｐｒａｌｌ、Ｄｉｇｉｔａｌ　Ｉ
ｍａｇｅ　Ｐｒｏｃｃｓｓｌｎｇ、３７ｇ−４２５，（
ｊｏｈｎ　　ＷｉＩａｙ　　＆　　５ｏｎｓ；１９７ｇ
）。

これに対して、コンボリューション法を採る場合は、装
置の構成のシンプルさや処理速度の速さといった利点を
生かせるように「縦３画素」　「横３画素Ｊ　　（３Ｘ
３）、「縦５画素」　「横５画素」（５Ｘ５）といった
小さなマスクサイズにより実現される。ところが、限ら
れたマスクサイズで文Ｍｌ）、２）に紹介されているよ
うなフィルタの特性を実現する方法については確立され
ていない。従来の２次元のコンボリューションフィルタ
（一般にはＦＩＲフィルターと呼ばれている）を設計す
る方法としては、下記文献３）に紹介されているマクレ
ラン変換か一般的である。

文献３）武部、　ディジタルフィルターの設計”５１−
９５．２１１−２１３．　　（東海大学出版会、１９８
６）。

これは１次元で設計したフィルターを２次元へ周波数変
換するものである。１次元フィルターの設計法としては
窓関数法１周波数サンプリング法など種々の方法が提案
されているが、中でも重み付きチエビシエフ近似法が代
表的である。これは、直線位相フィルターを仮定した場
合、伝達関数は余弦級数で表わされるが、その係数を、
近似誤差関数の絶対値の最大値が最小になるように求め
る方法である。この１次２元のチエビシエフ近似法と２
次元に変換するマクレラン変換とを組み合イ）せると、
かなり複雑なアルゴリズムとなる。特に設定するパラメ
ータが多くあり、その中にはかなり経験を積まないと旨
く設定できないものがある。

しかも低域通過フィルターなど単純な伝達特性を有する
ようなコンボリューションフィルターでないと適用が難
しい。

［発明が解決しようとする課題］ 従来の空間周波数フィルタリング装置のうち、目的実現
のための最適な処理を行なうためにはフーリエ変換法を
採用するのが望ましい。しかし装置を構成するにはＦＦ
Ｔ演算器や画像メモリが新たに必要で、規模が大きくな
ってしまう。その上フレームフレームで処理を完了させ
るようなリアルタイム処理は難しい。そこで装置の簡便
化や処理の高速化のためには、小さめのマスクサイズに
よるコンボリューション法を採用するのか望ましいとい
える。しかし従来の２次元コンボリューションフィルタ
ーの設計法は、アルゴリズムが複雑で、経験的な技術が
必要である上に、比較的単純な特性のフィルターでない
と適用が難しいという欠点があった。

そこで本発明は、空間周波数領域で設＝１したフィルタ
ーの特性に最適近似させたコンポリューンヨンフィルタ
ーを簡便な方法で設＝１し、コンボリューション法によ
り比較的小形で高速度な空間周波数フィルタリングを行
なえる装置を提供することを目的とする。

［課題を解決するための手段および作用］空間周波数フ
ィルタリング装置をコンホリュションフィルターの設計
部と、コンポリュージョン法によるフィルタリング部と
で構成する。そして設計部においては、まず空間周波数
領域でフィルターの設計を行ない、この特性に最適近似
したコンボリューションフィルターをモンテカルロ法に
より設、ｉｌする方法を提案する。

第１図はフィルター設＝１の手順を表わす流れ図を示す
。ここで扱うフィルターは等方向な特性をａする０位相
のものとする。第１図に示すように、まずフィルターを
空間周波数領域で設計し、これを理想フィルター／ｌ’
（、とする。Ｎｏを逆フーリエ変換し、Ｎｏのインパル
スレスポンスＮ。を求める。等方向な０位相のフィルタ
ーを仮定するとコンボリューションフィルターは、動径
方向の係数列ｋｕ　　（ｒ）だけを考慮すれば良い。Ｌ
ＸＬ　（Ｌは奇数）のサイズのコンボリューションフィ
ルターは（Ｌ−１）／２＋　１個の要素からなる係数列
ｋｋ＋（ｒ）により設定できる。本方法ではＭ　−（Ｌ
−１）／２＋２個の係数列ｋｖ　　（ｒ）を設定し、ｋ
ｖ　　（ｒ）−０とする。まず、最初に設定するＭｔ、
−（Ｌｔ、−１）／２＋２個の係数列ｋ　ｂｌ＋＋（ｒ
）を求める。表現を簡潔にするために、以下行列表現す
る。ｋＭ　（ｒ）をベクトルｌｋＭで表わす。

ＬｋＭ−［ｋ　（１）　　ｋ　（２）　　・・・ｋ（Ｍ
）］・・・・・・■ 初期設定の係数列ベクトルｋＭｏは次式で求める。

ｌｋＭｏ−（Ｓ　ＨＯＴＭＯ）“　・・・・・・・・・
・・・・・・■扱う画像の大きさ、っまりＨ８のサイズ
はＮＸＮとする。Ｓ、ＴＭＯは次式で表わされる。

Ｓ−［１，Ｏ，・・・、　　Ｃ１］・・・・・・■Ｉｋ
ＭＬｌのうち、Ｍ、１行目の係数はＯに置き換える。

次にＩｋ　ｙ＋ｏを基に２次元のコンボリューションフ
ィルターＨＬＯを設定する。

第２図はＨＬの設定の概念図である。原点を中心にＬＸ
Ｌの領域を設定し、各画素に設定する値は、原点からの
距離ｇに応じて、動径方向の係数列ｋＭ　（ｒ）から線
形補間により求める。こうして設定したコンボリューシ
ョンフィルターＨＬｏをフーリエ変換し、動径方向の空
間周波数り、。を算出する。

ｈＬｏ−（Ｓ　（ＦＨ＋、Ｆ）Ｔ　（Ｎ／２＋１））・
・・・・・■ ただしＦはフーリエ変換のオペレータ行列である。

１　１　ω０　　ω０　　・・・　　ω０Ｆ−□　１　
ωｔｌ　　　ω１　　　　　ωＮ−１１・・・■Ｊ““
““Ｉ’ｒｌω（ω２ω２（Ｎ−凰）ω（゛　　　ωＮ
−１ω（Ｎ−１＋ ただし、ω−ｅｘｐ［−２πｊ　／　Ｎ　］ＴＬＮ２゜
１．は０式と同様に定義される。

評価関数ｇ′は次式で定義される。

ｇ　　−１１ｈａ　−ｈＬ　Ｉｆ２−・−−−−■ただ
し、ｈｏは理想フィルター特性／／、の動径方向成分で
ある。

ｈｏ　＝　（Ｓ　ｉｎ　Ｔ（Ｎ　２−１１　）　’　・
・・・・・■最初の評価間数値は ｇ’　＝　Ｉｆ　ｈａ　　ｈＬｏｌｌ　２・・・・・・
［相］となり、ｇ−ｇ　　とおく。次にモンテカルロ法
により最適なコンボリューションフィルターの係数値を
求めていく。つまり、ある微小Ｈ（ｇｒａｉｎ　：）を
ｋＭｏの各要素値に加えて０式の評価関数値を旧算し、
ｇ＜ｇ’　となった場合はｇをｇ′に置き換え、ｇｒａ
ｊｎを受け入れる。ｇがある値に収束したと考えられる
時点て繰り返しを止め、最適なフィルターとする。こう
してり。ＸＬ、のサイズのコンボリューションフィルタ
ーの設計か終了すると、これから（Ｌ（、−２）Ｘ　（
Ｌｏ−２）の領域を切り取り、これを初期値として（Ｌ
Ｌ、　−２）　Ｘ　（Ｌ−２）のサイズのフィルターを
モンテカルロ法により求める。実際は、動径方向の係数
列ｌ（、。

（ｒ）の（Ｍ、、−１）番目の値を０に置き換え処理を
続行することになる。このようにしてコンボリューショ
ンフィルターのサイズを小さくしていき、最終的に所望
のサイズのフィルターを求める。

以上のような方法により、所望のサイズ、例えば５×５
なとの小さなサイズのコンボリューションフィルターを
求める場合に、より理想フィルター特性Ｎ。に近いフィ
ルターに比較的少ない繰り返し回数で近似させることが
できる。しかも手法は極めて簡便であり、紅験的技術や
特殊な手法を必要としない。

［実施例］ （第１実施例） 第３図は本発明の第１の実施例の全体構成を示す斜視図
である。不装置は図示の如くコンピュータ１００．コン
ボリューションフィルター設計装置２００９画像人力装
置３００．フィルタリング装置４００．ＴＶモニタ５０
０からなっている。

コンピュータ１００では、空間周波数帯域での理想フィ
ルター特性Ｎ（＋が算出される。その結果は、コンボリ
ューションフィルター設＝１装置２００に転送される。

コンボリューションフィルター設＝１装置２００では所
定のマスクサイズによる最適なコンボリューションフィ
ルターが導出される。そのコンボリューションフィルタ
ーの係数値は、フィルタリング装置４００に転送される
。フィルタリング装置４００では画像人力装置３Ｌ）０
により入力された画像信号に対し、コンボリューション
法による空間周波数のフィルタリングか実行される。そ
の処理画像はＴＶモニタ５００により表示される。本構
成においてコンボリューションフィルター設＝１装置２
００と、フィルタリング装置４００とは、オンラインで
接続されているか、コンボリューションフィルター設：
！゛装置２０（）により設計された係数値をフィルタリ
ング装置４０Ｕ内のＲＯＭに記録するぢしてオフライン
で接続される構成であっても良い。また、画像人力装置
３００は画像を人力すると同時に何等かの処理を行うよ
うな装置であっても良いし、フィルタリング装置４００
は空間周波数のフィルタリングを前処理として、その後
にさらに画像処理を行なって、結果をＴＶモニタ５００
に表示するようなＩｆ４成であっても良い。

第４図はコンボリューションフィルター設＝１装置２０
０の構成を示すブロック図である。ます、コンピュータ
１００により設計された理想フィルター７１１０は、イ
ンタフェース（１／Ｆ）２０．１とセレクタ２０２を介
して、メモリ２０３に一度記録される。理想フィルター
特性／／　（＋の動径ｈ゛向のプロフィールｈ。は、Ｎ
ｏの原点から横方向へＮ／２＋１個の数列で表現される
ので、この数列り。がメモリ２０３より読み出され、セ
レクタ２０４を介して、メモリ２１３に記録される。そ
して、メモリ２０３から読み出された理想フィルター特
性Ｎ。はセレクタ２０４を介してＦＦＴ演算器２０５に
人力され、フーリエ変換か実行され、結果Ｈｏはメモリ
２０３に再び記録される。■式で表わされる初期設定に
おける動径方向の特性＋１（ＭＯは、メモリ２０３に記
録されている理想フィルターのインパルスレスポンスＨ
６における原点から横方向へのＭ。個の数列で構成すれ
ば良い。したがってこの数列１ｆ（Ｍ（＋は、メモリ２
０３から読み出され、セレクタ２０４．セレクタ２０６
を介してメモリ２０７に記録される。モンテカルロ法に
よる近似工！算は次のような構成で行なう。メモリ２０
７に記録されている動径方向の特性１にいと、グレイン
発生器２０８により発生されたグレインとが、加Ｗ、／
減算器２０９により加算もしくは減算され、結果Ｄ（１
，はメモリ２１０に記録される。グレインを付加された
係数列１ｋＮ１’から２次元のコンボリューションフィ
ルターＨＬか先に述べた方法により設定される。これは
、メモリ２１０に記録されである係数列に１．′かセレ
クタ２］１を介して演算器２１２に必要に応じて入力さ
れ、演算器２１２において所定の補間演算が行なわれる
ことにより実現される。その結果は、セレクタ２０２を
介してメモリ２０３に記録される。メモリ２０３に記録
されたコンボリューションフィルターＨＬはＦＦＴ演算
器２０５によりフーリエ変換され空間周波数特性ＨＬが
算出される。その結果はメモリ２０３に再び記録される
。メモリ２０３からは、０周波数の値が読み出され、Ｃ
ＰＵ２２υに送られる。ＣＰＵ２２０には、理想フィル
ター／７’ｏの０周波数のフーリエ係数値Ｈ６（０）か
あらかしめ記録されており、新たな周波数特性／７′Ｌ
の０周波数のフーリエ係数値ＮＬ　（０）か／／、、（
０）と同じにｔよるような乗数Ｃか算出される。メモリ
２０３から読み出された動径方向のプロフィールめＬは
、セレクタ２０４を介して乗算器２１４に入力される。

そしてＣＰＵ２２０により算出された乗数Ｃと乗算され
た後に、減算器２１５においてメモリ２１３に記録され
ている理想特性ｌｈｏと減算か行なわれる。減算器２１
５からの出力Ｉｈ。

ｌｈｌ、は累積加算器２１６に人力され、成分の２乗和
、つまりノルムｇ　’　＝　ＩＩ　ｈ　ｕ　　ｌｈ　Ｌ
　ＩＩ　２が計算される。その結果はＣＰＵ２２０に転
送される。

ＣＰＵ２２０では記録されている評価関数値ｇと新たな
ｇ′との大小を調べ、もし、ｇ＞ｇ’であればメモリ２
１０に記録されている係数列１ｋ　Ｍをメモリ２０７に
記録されている係数列１ｃ　Ｍに置き換える。そして、
グレイン発生器２０８に新たなグレインを発生させ、上
記の処理を繰り返す。

この様にして、あるマスクサイズＬＸＬのコンボリュー
ションフィルターの近似計算が行なわれる。

そして評価関数値ｇの走化が十分小さくなった時点て繰
り返しを止める。なお解か収束すると思われる繰り返し
回数をあらかじめ設定しておき、所定の回数繰り返すよ
うにしても良い。ＬＸＬのマスクサイズでの処理か終わ
ると、次にサイズを（Ｌ−２）Ｘ　（Ｌ−２）に落とし
、新しいサイズの最適フィルターを求める。実際にはメ
モリ２０７に記録されている係数列１（ＭのＭ−１番目
の値を０に置き換え、モンテカルロ法による処理を続行
することになる。こうして、所定のマスクサイズに至る
まで処理を続け、最終的に求められたコンボリューショ
ンフィルターＨＬの各係数値かインタフェース２２１を
介してフィルタリング装置４００へ転送される。

第５図は、フィルタリング装置４００の構成を示すブロ
ック図である。コンボリューションフィルター設計装置
２００により設計されたフィルターの係数値はインタフ
ェース４０１を介してＲＯＭ４０２−１，４０２−２に
記録される。一方、ＴＶ左カメラの画像人力装ｆＥ３０
０により人力された画像信号はＡ／Ｄ変換器４０３によ
りディジタル信号に変換され、空間周波数フィルタリン
グが実行される。本実施例では、５×５のマスクサイズ
によるフィルタリングを仮定する。人力ディジタル画像
信号は、ラインバッファ４０４−１〜４０４−４を介し
て多重回路（ＭＵＸ）４０５−１．４０５−２に人力さ
れる。ＭＵＸ４０５−１に注目すると、図のように２つ
のラインバッファ４０４−１．４０４−２を用いること
により、同じ列の隣り合った３行の画像信号が同時に人
力されるように構成されている。ＭＵＸ４０５−１では
並列に入力される３行の画像信号を所定の逐次的な信号
列に並び換え、Ａ／Ｄ変換器４０３の速度の少なくとも
３倍以上の速度で出力する。例えばディジタル信号が６
ＭＨ２の速度で入力され、ＭＵＸ４０５−１は、２０Ｍ
Ｈ２で信号を出力するように設定される。ＭＵＸ４０５
−２はＭＵＸ４０５−１に入力される画像信号に隣接す
る下２行の画像信号を逐次系列に変換し、Ｍ　Ｕ　Ｘ　
４０５−１と同じ速度で出力する。Ｍ　Ｕ−Ｘ　４０５
−１　。

４０５−２からの出力信号はディジタルフィルタリング
プロセッサ（ＤＦＰ）４０６−１．４０６−２に人力さ
れる。ＤＦＰ４０６−１では、ＭＵＸ４０５−１からの
出力速度に同期してバイブライン方式により３行５列の
局所領域についてＲＡＭ４０２−１に記録されているフ
ィルタリング係数との積和演算が行なわれる。その結果
は人力信号の１／３の速度で出力され、レジスタ４０７
−１に保持される。ＤＦＰ４０６−２では同様に２×５
の局所領域の積和演算が行なわれ、結果はレジスタ４０
７−２に保持される。そして両レジスタ４０７−１．４
０７−２に保持されているデータは、加算器４０８によ
り加算され、Ａ／Ｄ変換器４０３と同じ速度のディジタ
ル信号として出力されてる。そしてＤ／Ａ変換器４０９
により所定のアナログビデオ信号に変換される。この様
にして、５×５の局所領域についてのコンボリューショ
ンフィルタリングが実行される。

上記した第１本実施例によれば、本発明の基本概念に基
づいてコンピュータ１００により設計された理想的フィ
ルター特性に最適近似させた所定のマスクサイズのコン
ボリューションフィルターを設計し、そのフィルターを
使って画像の空間周波数フィルタリングを行なう作用が
有る。特にフィルタリング装置４００では、ディジタル
回路によりリアルタイムで処理が行なわれる。かくして
フィルターの設計の自由度が大きく応用性の高い装置を
実現できる。またディジタル処理であるため再現性に優
れ、しかも高速に処理できる。

（第２実施例） 本実施例は本発明を画像処理による長焦点深度顕微鏡装
置へ応用したものである。本実施例は、解像度や明るさ
を損なうこと無く、顕微鏡画像の焦点深度を深くするこ
とを目的として、次のような画像処理手法を実行するも
のである。対象物の深さ方向の構造に応じて、観察すべ
き全ての物体面に合焦させるようにステージを駆動しな
がら画像を積算人力し、適当な中域〜高域強調フィルタ
リングを実行することにより、観察すべき全ての物体面
に同時にフォーカスの合った画像を再生する。なお、本
手法の概念は、先出願特開平０１−３０９４７８号に開
示されている。

本実施例の構成のうち、理想フィルター特性を算出する
ためのコンビニ−・夕と、コンボリューションフィルタ
ー設計装置とは、第１実施例に示したコンピュータ１０
０とコンボリューションフィルター設計装置２００と同
様である。

第６図は第２実施例の顕微鏡画像処理装置の構成を示す
図である。図に示すように本装置は顕微鏡６００．ＴＶ
カメラ７００．ブロセッ”ｊ　８００　。

ＴＶモニタ９００それにマンマシーンインタフェース１
０００とからなっている。顕微鏡装置６００のステージ
部には、Ｘ−Ｙステージ駆動装置６０１が設けられてお
り、プロセッサ８００からの指令によりステージのＸ−
Ｙ方向への移動がコントロールできるようになっている
。また、ステージ上には振動式ホルダー６０２が設けら
れており、圧電素子などの振動アクチュエータによりホ
ルダ上の対象物を、光軸方向に所定の振動幅でＴＶカメ
ラ７００のフレームもしくはフィールド蓄積時間よりも
十分速い同期で振動駆動する。プロセッサ８００１：　
ハ、Ａ／Ｄ変換器８０１　カあり、ＴＶカメラ７００に
より入力された画像信号をディジタル変換する。本装置
は本処理動作を行なう前に、処理条件の設定を行なうた
めの前処理を行なうものとなっている。

前処理時には、Ａ／Ｄ変換器８０１からのディジタル画
像信号はコントラスト検出部８２０に入力される。コン
トラスト検出部８２０は、図に示すようにバンドパスフ
ィルタ８２１１乗帥器８２２、加算器８２３、それにバ
ッファ８２４により構成されている。この検出部８２０
により人力画像のコントラストが検出され、これより振
動式ホルダ６０２の振動幅がＣＰＵ８５０において設定
される。つまり、振動式ホルダ６０２の振動幅を変えな
がら、入力画像のコントラストを検出し、コントラスト
値があるスレッショルド値以下になる振動幅を検知して
、これを本処理時における振動幅とする。この振動幅設
定法の詳細は前記先出願特許に記載されている。

本処理においては、所定の条件で光軸方向に振動された
対象物の画像信号がＴＶカメラ７００により入力される
と、この画像信号はプロセッサ８００内の、Ａ／Ｄ変換
器８０１によりディジタル信号に変換された後、画像処
理部８１０に入力される。Ｘ−Ｙステージ駆動装置６０
１によるステジ駆動が連続的に行なわれる場合は、人力
画像は画像処理部８１０においては、直接フレームメモ
リ８１２に記録される。そして引き続いて第１実施例の
フィルタリング装置４００と同様にラインバッファ　Ｍ
ＵＸ、ＤＦＰ、　　レジスタ、加算器で構成されるフィ
ルタリング装置８１３に人力されて、所定の空間周波数
フィルタリングが行なわれる。なお、コンボリューショ
ンフィルターの係数値はあらかじめコンピュータ１００
およびコンボリューションフィルター設計装置により設
定されており、ＲＯＭ８１４に記録されている。そして
、前処理により設定された最適振動幅に基づいて、ＲＯ
Ｍ８１４に記録されている係数のうち適当な係数が読み
出される。また、Ｘ−Ｙステージを静止させて観察する
場合は、加算２７１８１１とフレームメモリ８１２によ
り数枚から数１０枚の入力画像を累積加算し、−層Ｓ／
Ｎに優れた画像にしてから、フィルタリング処理を行な
う。このようにして人力され、フィルタリング処理を施
された画像は、Ｄ／Ａ変換器８０２によりアナログビデ
オ信号に変換されてＴＶモニタ９００に表示される。な
お、Ｘ−Ｙステージ駆動装置６０１はＸＹステージ駆動
ドライバ８４０により駆動され。

振動式ホルダ６０２は振動アクチュエータドライバ８４
０により駆動制御される。両ドライバはＣＰＵ８５０に
より制御される。また条件の設定等は、ＣＰＵ８５０に
接続されたマンマシーンインタフェース１０００により
観察者が行なうように構成されている。

理想フィルター設＝１のためのコンピュータ１００では
、積算入力した画像から全ての物体面にフォーカスの合
った画像を画成する回復処理のフィルター特性が以下の
ように算出される。ます、使用する対物レンズの特性か
らある物体面と合焦面との距離（Ｚ）より、ぼけによる
レスポンス関数（ＯＴＦ）が計算される。そしてこの計
算結果からある物体面の範囲Ｚに渡って積や人力した画
像のほけ息、つまり積算人力のレスポンス関数〃（ｕ、
ｖ）がｐ出される。なお（ｕ、ｖ）は空間周波数面の直
交座標である。

Ｄ（ｕ、ｖ）−ｆｚ　ｐ（ｕ、ｖ、ｚ）ｄｚ　　　　　
＝−＠たたし、Ｄ（ｕ、ｖ；ｚ）は、焦点はずれの距離
か２の場合のＯＴＦである。

・・・０ また、（ｘ、ｙ）は２次元の空間座標、Ｐ　（ｘ。

ｙ、ｚ）は焦点はずれ２の場合の一般化瞳関数であり、
０式で表わされる。λは設定した光の波長、ｆはレンズ
の焦点距離である。

Ｐ（ｘ、ｙ；ｚ）−Ｐ（ｘ、ｙ）ＥＸＰ［ｊｋＷ（旧Ｚ
）］　　　−＠ただし、Ｐｏ　（ｘ、ｙ）は合焦時の瞳
関数で０式で表わされる。

１−　　Ｘ　　＋ｙ　　　　　　　　　　−Ｏｋは波数
であり、ｋとλは、次式の関係がある。

２π λ Ｗ（ｒ７ｚ）は波面収差で、高倍率の対物レンズのよう
にＮ、Ａが大きいレンズの場合、０式で近似される。

ｒ”ｚ Ｗ　（ｒ　；ｚ　）−・・・０ ２　　（ｆ　２　＋ｒ　２　） 以上のような関係により積算画像のＯＴ　Ｆ　ｌ）　（
ｕ　。

Ｖ）がコ１°算される。一方、合焦時のＯＴＦ／、。

（ｕ、ｖ）は０式で算出される。

・・・０ したがって、回復処理のフィルター特性／／　（Ｕ　。

■）は最も単純には次式で表わされる。

１０　（ｕ、ｖ） Ｈ（ｕ、ｖ）−・・・Ｏ ｐ（υ、ｖ） 本実施例を半導体ｔＣの検査装置へ適用する場合のよう
−に、対象物が限定され、画像およびノイズの統計的性
質が予測＝Ｊ能な場合、０式で表わされるようなウィー
ナフィルターによりノイズを低減しながら回復処理を行
なわせるようにしても良い。

／／ｗ　（ｕ、ｖ）＝［ム、（ｕ、ｖ）　　　ｉ　ｐ（
ｕ、ｖ）　　ｌ　２］／　［ｐ（ｕ、ｖ）１１　　ｐ（
ｕ、ｖ）　　ｌ　２＋　Ｓ　　−−（ｕ、ｖ）／Ｓ　　
ｇｇ（ｕ、ｖ）１１・・・０ ただし、５−（Ｌｌ、　　Ｖ）はノイズの、Ｓ、、（ｕ
。

■）は理想画像の統計的に予測されるノくワースベクト
ルである。

このようにして設計された理想フィルター特性／／（Ｕ
、Ｖ）に近似させたコンボリューションフィルターをコ
ンボリューションフィルター設計装置により導出する。

なおコンボリューションフィルターは積算人力の範囲２
として想定される数段階の値を仮定して算出し、プロセ
ッサ８００内のＲＯＭ８１４に記録する。そして実際の
処理では、前処理によって設定した積算範囲に最も近い
条件のフィルターを読み出すことになる。

こうして設計されたコンボリューションフィルターによ
り、前記概念に基づく画像処理が実現される。

かくして第２実施例によれば、前記先出願特許に基づく
長焦点深度顧、微鏡装置を、より実用的な構成を有する
装置となし得る。特に空間周波数領域で設計した回復処
理のフィルターに最適近似したコンボリューションフィ
ルターを簡単に算出でき、最も好ましい処理画像を簡便
な構成で迅速に表示することが可能になる。

なお顕微鏡以外にも電子カメラや内視鏡など種々の光学
機器について、同様に応用することが可能である。

（第３実施例） 本実施例は、モンテカルロ法をＸＩ　Ｈする際に、フィ
ルターの特性だけを考慮するのではなく、推定観測画像
と実際の観測画像との２乗誤差を評価基準として最適な
コンボリューションフィルターを設計したものである。

第７図は全体構成を示す図である。図示の如く本装置は
入力系の劣化関数や理想特性を記憶し、フィルター特性
の初期設定を算出するコンピュータ１１００と、コンポ
リューシタンフィルター設計部１２００．画像入力装置
ｔ１３００．画像人力装置１３００からの画像信号を分
配する分配器１４００、フィルタリング装置１５００お
よびＴＶモニタ１６００から構成される。フィルタリン
グ装置１５００は、第１実施例に示したフィルタリング
装置＆４００と同様に構成されている。

第８図はコンポリューシジンフィルター設計部１２００
の構成を示す図である。コンピュータ１１００で設定さ
れたフィルター特性の初期設定ＩＯは、インタフェース
１２０１．セレクタ１２０２を介して、メモリ１２０３
に記録される。メモリ１２０３に記録されたフィルター
特性ＩｆｏはＦＦＴ演算器１２０６により２次元フーリ
エ変換が実行される。そしてインパルスレスポンス（点
像分布関数：ＰＳＦ）の初期設定Ｈ６が算出されて、メ
モリ１２０３に再び記録される。またＭ　ｌ）Ｊのマス
クサイズＬ。ｘＬｏのコンボリューションフィルターを
設定するための、Ｍ　（＋　”　Ｌ　（１／　２　＋２
個の要素からなる動径方向の係数列ｋｖａ（ｒ）か読み
出される。これはセレクタ１２０５、セレクタ１２０７
を介してメモリ１２０８に記録される。

次にコンピュータ１１００からは画像入力装置１３００
の劣化関数〃が読み出される。そしてインタフェース１
２０１およびセレクタ１２ｏ２を介して、メモリ１２０
４に記録される。また、画像入力装置１３００により入
力されたサンプル画像■がＡ／Ｄ変換器１２１４により
ディジタル画像信号に変換される。そしてセレクタ１２
ｏ２を介して、メモリ１２０３に記録される。同時に画
像信号Ｖ　ハ、セレクタ１２ｏ５を介してメモリ１２１
５に記録される。そしてメモリ１２０Ｂに記録された画
像信号ＶはＦＦＴ演算器１２０６により２次元フーリエ
変換が行われる。その結果は再びメモリ１２０３に記録
される。メモリ１２０３に記録されているサンプル画像
のフーリエ変換γはセレクタ１２０５を介して、またメ
モリ１２０４に記録されている劣化関数ｐはセレクタ１
２１６を介して、乗算器１２１７に入力される。そして
両者は各空間周波数成分ごとに乗算されて、結果〃・γ
はメモリ１２０４に記録される。このようにして初期設
定が完了し、次にモンテカルロ法による近似計算か以下
のように行なわれる。メモリ１２０８に記録されている
動径方向の係数列ｋ　ｈ＋と、グレイン発生器１２０９
により発生されたグレインとが加算／減算器１２１０に
より加算もしくは減算される。その結果に２．′はメモ
リ１２１１に記録される。メモリ１２１１に記録された
係数列ｋＭ　は捕間演算器１２１３により２次元のコン
ボリューションフィルターＨＬに変換される。

その結果はセレクタ１２０２を介してメモリ１２０３に
記録される。そして、コンボリューションフィルターＨ
ＬはＦＦＴ演算器１２０６により２次元フーリエ変換さ
れる。その結果は再びメモリ１２０３に記録される。こ
のメモリ１２０３に記録されたコンボリューションフィ
ルターＨＬの空間周波数７／Ｌのうち０周波数の値か読
み出され、ＣＰＵＩ　２３０に入力される。ＣＰＵＩ　
２３０では、空間周波数Ｉの０周波数値が常に同じ値に
なるような係数Ｃを算出し、この値をセレクタ１２１６
を介して乗算器１２１７へ送る。乗算器１２１７ては、
メモリ１２０３から読み出された空間周波数特性ＮＬの
各空間周波数成分値と係数値Ｃとが乗算される。その結
果ＩＬ　　はメモリ１２０３に記録される。次にメモリ
１２０３から空間周波数特性ＨＬ、メモリ１２０４から
先に計算されたｐ・γが読み出される。そして乗算器１
２１７によって両者が空間周波数ごとに乗算され、結果
ｐ・γ・ＮＬ′はメモリ１２０３に記録される。

このメモリ１２０３に記録された空間周波数特性はＦＦ
Ｔ演算器１２０６により２次元逆フーリエ変換される。

そして推定観測画像Ｖがメモリ１２０３に記録される。

次に推定観測画像Ｖは、セレクタ１２０５を通過した後
にセレクタ１２１８を介してレジスタ１２１９−１に保
持されると同時にレジスタ１２１９−２に保持される。

レジスタ１２１９−１および１２１９−２に保持された
値は、乗算器１２２０により乗算されてから、加算器１
２２１およびバッファ１２２２により累積加算される。

この様にして、推定観測画像■の内積Ｖ′Ｖが計算され
、その結果はセレクタ１２２３を介してレジスタ１２２
４に記録される。次にメモリ１２０３に記録されている
推定観測画像Ｖと、メモリ１２１５に記録されている観
測画像Ｖとは乗算器１２２０により乗算される。そして
加算器１２２１およびバッファ１２２２により累積加算
される。この様にして推定観測画像Ｖと観測画像Ｖとの
内積Ｖｌｖが計算され、その結果はセレクタ１２２３を
介して乗算器１２２５に入力される。

そしてＣＰＬ１１２３０から読み出された係数値゛２゛
と乗ｐされた後、減算器１２２６においてレジスタ１２
２４に記録されている推定観１ｌｐ１画像■の内積■１
ｖと減算される。この禄にして評価関数値ｇ′−”Ｆｏ
Ｖ−２ｖ“Ｖが計算され、その結果はＣＰＵ１２３０に
入力される。ＣＰＵ１２３０では、記録されている評価
関数値ｇと新たに人力された値ｇ′との大小関係を調べ
、もしｇ〉ｇ　であればメモリ１２１〕に記録されてい
る係数列ｋＭ　をメモリ１２０８に記録されている係数
列ｋＭに置き換える。そして、グレイン発生器１’２０
９に新たなグレインを発生させ、上記の処理を繰り返す
。以上のような構成により初期設定のり。ｘ　Ｌ　ｏの
マスクサイズのコンボリューションフィルターからサイ
ズを小さくしていき、最終的に所望のサイズのコンポリ
ューシぢンフィルターを決定するのは、ｔｔ５１の実施
例に示した通りである。最終的に求められたコンボリュ
ーションフィルターの各係数値はインタフェース１２２
７を介して、フィルタリング装置１５００へ転送される
。

第１．２実施例で示した設計法では評価関数を［相］式
で定義したが、これは単にコンボリューションフィルタ
ーの空間周波数特性が理想特性に近くなるように２乗誤
差を定義したもので各空間周波数における誤差が実際の
処理画像にどう影響するかまでは考慮されていない。そ
こで、劣化画像の復元のような目的にフィルタリングを
行なう場合に、評価関数を観測画像と推定観測画像との
２乗：！′１差で定義する。推定観測画像とは、観測画
像をコンボリューションフィルターにより復元すること
により原画像を推定し、さらにこれを劣化関数で劣化さ
せたものをいう。この推定観測画像と実際の観測画像と
の２乗誤差が小さくなるようにモンテカルロ法を適用す
ることにより、最適な推定）Ｍ画像を復元するようなコ
ンボリューションフィルターを導出することができる。

本実施例の方法は、画像の種類がある程度限定され、劣
化関数ｐおよび劣化の少ない理想伝達関数ノ、１が分か
っている場合に後述するような統計的手法によりコンボ
リューションフィルターを求めるものである。

これは例えば第２実施例で示した長焦点深度顕微鏡装置
で、半導体ＩＣを観察するような場合に有効である。そ
の場合、劣化関数りは、合焦面を変えながら入力した画
像を加え合わせる操作によるＯＴＦ、理想伝達関数（理
想的０ＴＦ）／、、は、用いた対物レンズで合焦させた
際のＯＴＦである。

観測面°像Ｖはサンプルとして用いた半導体ＩＣの入力
加算画像を用いることになる。

以下に本実施例の作用を詳述する。２乗誤差ｅ２は次式
で定義される。

ｅ２−Ｅ　　（［ｖ−ｖＦ”］　　’　　［Ｖ−’Ｖｌ
　１ｍｌ：　（Ｖ’　Ｖ）−２Ｅ　ＩＶ’　Ｖｌ十Ｅ（
ｖｌ？）・・・０ ただし、Ｖ：観測画像 Ｖ：推定観測画像 Ｅ（・）　　：％合平均のオペレータ であり、■およびＶは画素数をｎ２　（＝ｎＸｎ）とす
ると、ｎ　”　ｘ　１行列つまり要素数ｎ２の列ベクト
ルで表わされる。また推定観測画像Ｖは次式％式％ ただし、Ｄ−Ｆ”　ＤＦ”　：劣化関数ｐのＰＳＦＨＬ
−Ｆ’　／／ＬＦ”　：コンボリューションフィルター
　ＮＬはフンボ リューションフィルターＨＬ の空間周波数特性 γ−ＦＶＦ　：観測画像γのフーリエ変換ＦｉｌＬｅｒ
（Ａ−８）：空間周波数領域でＡとＢとの各空間周波数
の係数値を 掛ける。つまりフィルタリン グのオペレータ である。０式により０式の第１項は、コンボリューショ
ンフィルターＨＬに無関係なので、評価関数ｇ′は次式
のように残り２項で定義される。

ｇ’　−Ｅ　（Ｖ”　？ｊ　−２Ｅ　（Ｖ’　Ｖｔ−＠
０式の評価関数を用いることにより、より確からしい原
画像を推定するコンボリューションフィルターを設計す
ることができる。所望のマスクサイズのフィルターを設
計するまでの過程は、本発明の基本概念に基づくもので
あるが、初期値として□は、次式で表わされるＨ６のイ
ンパルスレスポンスＨ６の初期設定のマスクサイズＬ（
１×Ｌｏにおける動径方向の係数列ＩＩ（ＭＬＩ　　Ｍ
　（、つＬ　ｔ、／　２　＋２を用いる。

］（１ ｋ　Ｍｌ＋−（Ｓ−Ｈｕ　　・Ｔ　ＭＵ）＝　（Ｓ−Ｆ
”　・Ｎｏ　−Ｆ”　・ＴＭｏ）　−Ｑ本実施例の概念
は、画像の統計的性質を用いて、最適なフィルター設計
を行なうものであるので、複数のサンプル画像を用いて
０式を計算するのが望ましいが、サンプル画像は１枚で
あっても良い。

以上の説明では、１枚のサンプル画像からコンボリュー
ションフィルターを設計する場合について述べた。もち
ろん複数のサンプル画像を用いても良く、その場合は、
第８図にボしたコンボリューションフィルター設計部１
２００にｖＩｖおよびＶ゛■を異なるサンプル画像につ
いて計算し、それらを累積加算できるように構成を追加
するだけで実現できる。第９図は上記第３実施例の構成
に基づくフィルタ設計の流れを示す図である。

第３実施例によれば、扱う画像の性質を考慮することに
より、劣化画像の復元に最適なコンボリューションフィ
ルターを設計することかでき、コンボリューション法に
よる簡便なフィルタリング装置を用いて、最も確からし
い原画像を高速に復元することができる。特に、あらか
しめ理想フィルター特性を求めておく必要がなく、コン
ボリューションフィルターの設計の際に、フィルター特
性の最適化をも同時に実行できるという特徴を有する。

［発明の効果］ 本発明によれば、所望のマスクサイズで目的を達成する
のに最適なコンボリューションフィルターを簡便な方法
で設計することか可能であり、コンボリューション法に
より比較的小規模な装置て高速度に空間周波数フィルタ
リングを行なうことのできる空間周波数フィルタリング
装置を提供することができる。

【図面の簡単な説明】

第１図は本発明のフィルター装置の設計手順を示す流れ
図、第２図は２次元コンボリューションフィルターＨＬ
の設定の概念図である。第３図は本発明の第１実施例の
全体構成を示すブロック図、第４図は第３図のコンボリ
ューション設計装置の構成を示す図、第５図は第３図の
フィルタリング装置の構成を示す図、第６図は本発明の
第２実施例の全体構成を示す図、第７図は本発明の第３
実施例の全体構成を示す図、第８図は第７図のコンボリ
ューションフィルター設計部の構成を示す図、第９図は
第３実施例に係る複数サンプル画像フィルター設計手順
を表わす流れ図である。 １００・・・コンピュータ、２００・・・コンボリュー
ションフィルター設計装置、３００・・・画像入力装置
、４００・・・フィルタリング装置、５００・・・ＴＶ
モニタ。 第１図

Claims (3)

【特許請求の範囲】
1. （１）画像の空間周波数のフィルタリングを行なう装置
   において、 フィルターの理想的空間周波数特性を設定する手段と、 前記フィルターの実空間における動径方向の係数列を所
   定の長さで打ち切る手段と、 前記フィルターの実空間における動径方向の係数列に微
   小な数量を加える手段と、 前記フィルターの実空間における動径方向の係数列を実
   空間における２次元の係数行列に変換する手段と、 前記２次元の係数行列の空間周波数特性が、前記理想的
   空間周波数特性に近いか否かの評価基準を計算する手段
   と、 この手段により計算された前記評価基準から前記フィル
   ターの実空間における動径方向の係数列に加えられた微
   小な数量を受け入れるか否かを決定し、前記微小な数量
   を受け入れるかもしくは除外する手段と、 前記フィルターの実空間における動径方向の係数列に微
   小な数量を加えてから、その微小な数量を受け入れるか
   もしくは除外するまでの過程を繰り返す手段と、 前記フィルターの実空間における動径方向の係数列を所
   定の長さで打ち切る手段により係数列の長さを小さくし
   ながら、前記フィルターの実空間における動径方向の係
   数列に微小な数量を加えてから、その微小な数量を受け
   入れるかもしくは除外するまでを繰り返す過程を繰り返
   す手段と、前記全手段により最終的に求めた実空間にお
   けるフィルターの２次元の係数行列を用いて、画像の各
   画素の近傍領域において畳み込み演算を行なう手段と、 を有することを特徴とする画像の空間周波数フィルタリ
   ング装置。
2. （２）前記フィルターの理想的空間周波数特性を設定す
   る手段は、 前記フィルターの実空間における動径方向の係数列に微
   小な数量を加えてから、その微小な数量を受け入れるか
   もしくは除外するまでの過程を繰り返す手段により実現
   されるものであることを特徴とする請求項（１）に記載
   の空間周波数フィルタリング装置。
3. （３）前記２次元の係数行列の空間周波数が、前記理想
   的空間周波数特性に近いか否かの評価基準を計算する手
   段は、 前記２次元の係数行列の空間周波数特性が、サンプル画
   像に対して最適な処理画像が得られる前記理想的空間周
   波数特性に近いか否かの評価基準を計算する手段である
   ことを特徴とする請求項（１）に記載の空間周波数フィ
   ルタリング装置。