周波数領域での線形フィルターの設計

このトピックでは、周波数領域でフィルター処理を実行する関数について説明します。空間領域でのフィルターの設計方法の詳細については、空間領域でのイメージのフィルター処理とはを参照してください。

2 次元有限インパルス応答 (FIR) フィルター

Image Processing Toolbox™ ソフトウェアは、1 つのクラスの線形フィルターをサポートしています。それは、 2 次元有限インパルス応答 (FIR) フィルターです。FIR フィルターは、単点またはインパルスに対して有限の範囲をもっています。Image Processing Toolbox のすべてのフィルター設計関数は FIR フィルターを返します。

FIR フィルターは、MATLAB^® 環境で最適なイメージ処理を行ういくつかの特徴をもっています。

FIR フィルターは、係数の行列として表現することが簡単です。
2 次元の FIR フィルターは、1 次元の FIR フィルターを自然に拡張したものです。
FIR フィルター設計には、よく知られた信頼性の高い方法が複数存在します。
FIR フィルターは実装が簡単です。
FIR フィルターは、歪みを避けるために、線形位相で設計できます。

別のクラスのフィルター、無限インパルス応答 (IIR) フィルターは、それほどイメージ処理アプリケーションには適していません。FIR フィルターの設計や実現が簡単で、安定性があるのに対し、IIR フィルターは、この部分に欠陥があります。そのためこのツールボックスでは、IIR フィルターをサポートしていません。

メモ

この節で説明しているほとんどの設計法は、Signal Processing Toolbox™ の関数で作成される 1 次元のフィルターまたはウィンドウから 2 次元フィルターを作成する方法を使用しています。このツールボックスは必須ではありませんが、Signal Processing Toolbox ソフトウェアがない場合はフィルターの設計が困難になることがあります。

1 次元フィルターの周波数変換を使用した 2 次元フィルターの作成

この例では次を使用します。

ライブスクリプトを開く

この例では、関数ftrans2を使用して 1 次元 FIR フィルターを 2 次元 FIR フィルターに変換する方法を説明します。特定の特性をもつ 1 次元フィルターは、対応する 2 次元フィルターよりも簡単に設計できるため、この関数は役に立ちます。周波数変換法では、特に遷移帯域幅とリップル特性など、1 次元フィルターのほとんどの特性が保存されます。1 次元周波数応答の形状は、2 次元応答でも明らかです。

この関数は、周波数変換を定義する一連の要素である "変換行列" を使用します。この関数の既定の変換行列は、ほぼ円対称なフィルターを作成します。自分で変換行列を定義することで、別の対称を得ることが可能です (詳細については、Jae S. Lim,『Two-Dimensional Signal and Image Processing』, 1990 を参照してください)。

Signal Processing Toolbox の関数 firpm を使用して 1 次元 FIR フィルターを作成します。

b = firpm(10,[0 0.4 0.6 1],[1 1 0 0])

b = 1×11

    0.0537   -0.0000   -0.0916   -0.0001    0.3131    0.4999    0.3131   -0.0001   -0.0916   -0.0000    0.0537

1 次元フィルターを 2 次元フィルターに変換します。

h = ftrans2(b);

フィルターの周波数応答を表示します。

[H,w] = freqz(b,1,64,"whole");
plot(w/pi-1,fftshift(abs(H)))
title("One-Dimensional Frequency Response")

freqz2(h,[32 32])
title("Corresponding Two-Dimensional Frequency Response")

周波数サンプリング法を使用したフィルターの作成

ライブスクリプトを開く

この例では、周波数サンプリング法を使用して、必要な周波数応答をベースに 2 次元フィルターを作成する方法を説明します。

周波数サンプリング法では、周波数応答の形状を設定する点を行列として与え、これらの点を周波数応答が通るようにフィルターを作成します。周波数サンプリングでは、与えられた点の間での周波数応答の動作に制約がありません。通常、応答はこの部分でリップルになります。リップルは定数値周りの振動量です。実際のフィルターの周波数応答は、理想的なフィルターの周波数応答がフラットになる部分で、リップルをもつことがよくあります。

ターゲットとする、11 行 11 列フィルターの 2 次元周波数応答を定義し、表示します。

Hd = zeros(11,11); 
Hd(4:8,4:8) = 1;
[f1,f2] = freqspace(11,"meshgrid");
mesh(f1,f2,Hd)
title("Desired Frequency Response")

関数 fsamp2 を使用して、ターゲットの周波数応答に基づいてフィルターを作成します。fsamp2 は、入力行列 Hd の点を通過する周波数応答をもつフィルター h を返します。

h = fsamp2(Hd);

フィルターの周波数応答をプロットします。

freqz2(h,[32 32])
title("Actual Frequency Response")

実際の周波数応答のリップルを理想の周波数応答と比較して見てください。これらのリップルは、周波数サンプリング設計法の基本的な問題です。これらのリップルは、必要な応答内の鋭い遷移部で必ず生じます。

フィルターを大きくすればするほど、リップルの空間的な広がりは小さくなります。しかしフィルターを大きくすると、リップルの高さは低くなりますが、フィルター処理の計算時間は多くなります。必要な周波数応答に対してよりスムーズな近似を行うには、周波数変換法またはウィンドウ処理法を使用することを検討してください。

ウィンドウ法を使用したフィルターの作成

この例では次を使用します。

ライブスクリプトを開く

ウィンドウ処理法は、理想的なインパルス応答にウィンドウ関数を乗算して対応するフィルターを生成し、理想的なインパルス応答にテーパーを適用します。周波数サンプリング法と同様、ウィンドウ処理法は必要な周波数応答を近似するフィルターの周波数応答を作成します。しかしウィンドウ処理法は、周波数サンプリング法よりも良い結果を出す傾向があります。

このツールボックスには、ウィンドウベースのフィルター設計を使用して 2 次元フィルターを設計する 2 つの関数が用意されています。

fwind1は、指定された 1 つまたは 2 つの 1 次元ウィンドウを 2 次元ウィンドウに変換します。
fwind2は、指定された 2 次元ウィンドウを使用します。

fwind1 は、2 次元ウィンドウを作成するための 2 つの異なる方法をサポートしています。

単一の 1 次元ウィンドウを指定すると、fwind1 は回転と同様の処理を使用して 1 次元ウィンドウを変換します。2 次元ウィンドウはほぼ環状に対称になります。
2 つの 1 次元ウィンドウを指定すると、fwind1 はそれらの外積を計算します。2 次元ウィンドウは四角形かつ分離可能になります。

この例では、関数 fwind1 を使用して、単一の 1 次元ウィンドウから 2 次元ウィンドウを作成します。まず、11 行 11 列フィルターのターゲット 2 次元周波数応答を定義して表示します。

Hd = zeros(11,11);
Hd(4:8,4:8) = 1;
[f1,f2] = freqspace(11,"meshgrid");
mesh(f1,f2,Hd)
axis([-1 1 -1 1 0 1.2])
title("Desired Frequency Response")

1 次元ウィンドウを指定します。この例では、Signal Processing Toolbox の関数 hamming を使用します。

win1d = hamming(11);

目的の周波数応答 Hd と 1 次元ウィンドウから 11 行 11 列のフィルターを作成します。

h = fwind1(Hd,win1d);

フィルターの周波数応答を表示します。

freqz2(h,[32 32])
axis([-1 1 -1 1 0 1.2])
title("Actual Frequency Response")

参考