最新のリリースでは、このページがまだ翻訳されていません。 このページの最新版は英語でご覧になれます。

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

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

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 の関数で作られるウィンドウで 2 次元フィルターを作成する方法を使用しています。このツールボックスは必須ではありませんが、Signal Processing Toolbox ソフトウェアがない場合、フィルターを設計することは容易ではありません。

1 次元 FIR フィルターの 2 次元 FIR フィルターへの変換

この例では、関数 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 =

  Columns 1 through 9

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

  Columns 10 through 11

   -0.0000    0.0537

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

h = ftrans2(b);
h =

  Columns 1 through 9

    0.0001    0.0005    0.0024    0.0063    0.0110    0.0132    0.0110    0.0063    0.0024
    0.0005    0.0031    0.0068    0.0042   -0.0074   -0.0147   -0.0074    0.0042    0.0068
    0.0024    0.0068   -0.0001   -0.0191   -0.0251   -0.0213   -0.0251   -0.0191   -0.0001
    0.0063    0.0042   -0.0191   -0.0172    0.0128    0.0259    0.0128   -0.0172   -0.0191
    0.0110   -0.0074   -0.0251    0.0128    0.0924    0.1457    0.0924    0.0128   -0.0251
    0.0132   -0.0147   -0.0213    0.0259    0.1457    0.2021    0.1457    0.0259   -0.0213
    0.0110   -0.0074   -0.0251    0.0128    0.0924    0.1457    0.0924    0.0128   -0.0251
    0.0063    0.0042   -0.0191   -0.0172    0.0128    0.0259    0.0128   -0.0172   -0.0191
    0.0024    0.0068   -0.0001   -0.0191   -0.0251   -0.0213   -0.0251   -0.0191   -0.0001
    0.0005    0.0031    0.0068    0.0042   -0.0074   -0.0147   -0.0074    0.0042    0.0068
    0.0001    0.0005    0.0024    0.0063    0.0110    0.0132    0.0110    0.0063    0.0024

  Columns 10 through 11

    0.0005    0.0001
    0.0031    0.0005
    0.0068    0.0024
    0.0042    0.0063
   -0.0074    0.0110
   -0.0147    0.0132
   -0.0074    0.0110
    0.0042    0.0063
    0.0068    0.0024
    0.0031    0.0005
    0.0005    0.0001

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

[H,w] = freqz(b,1,64,'whole');
colormap(jet(64))
plot(w/pi-1,fftshift(abs(H)))
figure, freqz2(h,[32 32])

1 次元の周波数応答

対応する 2 次元の周波数応答

周波数サンプリング法

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

ツールボックス関数 fsamp2 は、2 次元の FIR フィルターに対して周波数サンプリング設計を実行します。fsamp2 は、入力行列 Hd の点を通過する周波数応答をもつフィルター h を返します。下の例では fsamp2 を使用して 11 行 11 列のフィルターを作成し、結果として作成されたフィルターの周波数応答をプロットします (この例の関数 freqz2 は、フィルターの 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]), colormap(jet(64))
h = fsamp2(Hd);
figure, freqz2(h,[32 32]), axis([-1 1 -1 1 0 1.2])

必要な 2 次元周波数応答 (左) と実際の 2 次元周波数応答 (右)

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

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

ウィンドウ法

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

ツールボックスには、ウィンドウをベースにしたフィルター設計に 2 つの関数 fwind1fwind2 が用意されています。fwind1 は、ユーザーが指定した 1 次元または 2 次元のウィンドウを使うことにより、2 次元フィルターを設計します。fwind2 は、指定した 2 次元ウィンドウを直接使用して 2 次元フィルターを設計します。

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

  • 1 つの 1 次元ウィンドウを変換し、回転と同様の処理を使用して、ほぼ環状に対称な 2 次元ウィンドウを作成します。

  • ウィンドウを 2 つの 1 次元ウィンドウに分離して、外積を計算することにより四角形を作成します。

次の例は、fwind1 を使用して、必要な周波数応答 Hd から 11 行 11 列のフィルターを作成します。この例では Signal Processing Toolbox の関数 hamming を使用して、1 次元ウィンドウを作成します。fwind1 はこれを 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]), colormap(jet(64))
h = fwind1(Hd,hamming(11));
figure, freqz2(h,[32 32]), axis([-1 1 -1 1 0 1.2])

必要な 2 次元周波数応答 (左) と実際の 2 次元周波数応答 (右)

必要な周波数応答行列の作成

フィルター設計関数 fsamp2fwind1、および fwind2 は、必要な周波数応答の大きさの行列をベースにフィルターを作成します。周波数応答は、さまざまな入力周波数に対するフィルターのゲインを記述する数学関数です。

必要な周波数応答行列を作成するには、関数 freqspace を使用します。freqspace は、任意の大きさの応答に対して等間隔の正しい周波数値を返します。freqspace によって返されるもの以外の周波数点を使用して必要な周波数応答行列を作成すると、非線形位相のような予期せぬ結果が生じます。

たとえば、0.5 のカットオフを持つ円の理想的なローパス周波数応答を作成するには、以下を使用します。

[f1,f2] = freqspace(25,'meshgrid');
Hd = zeros(25,25); d = sqrt(f1.^2 + f2.^2) < 0.5;
Hd(d) = 1;
mesh(f1,f2,Hd)

理想的な円のローパス周波数応答

周波数応答に対して、fsamp2fwind1fwind2 によって作成されるフィルターは実数になることに注意してください。この結果は、ほとんどの画像処理アプリケーションに適しています。一般的にこれを達成するために必要な周波数応答は、周波数の原点 (f1 = 0, f2 = 0) に関して対称になります。

フィルターの周波数応答の計算

関数 freqz2 は、2 次元フィルターに対する周波数応答を計算します。出力引数を設定しないと、freqz2 は、周波数応答のメッシュ プロットを作成します。たとえば、次の FIR フィルターについて考えてみましょう。

h =[0.1667    0.6667    0.1667
    0.6667   -3.3333    0.6667
    0.1667    0.6667    0.1667];

次のコマンドは、h の 64 x 64 点の周波数応答を計算し、プロットします。

freqz2(h)

2 次元フィルターの周波数応答

周波数応答行列 H と周波数点ベクトル f1 および f2 を得るには、出力引数を使用します。

[H,f1,f2] = freqz2(h);

freqz2 は、1.0 の値がサンプリング周波数の 1/2、あるいは π ラジアンに対応するように周波数 f1f2 を正規化します。

簡単な mn 列の応答に対しては、上に示すように、freqz2 は 2 次元の高速フーリエ変換関数 fft2 を使用します。また任意の周波数点ベクトルも指定できますが、この場合、freqz2 は時間がかかるアルゴリズムを使用します。

高速フーリエ変換、および線形フィルター処理、フィルター設計への応用の詳細は、フーリエ変換を参照してください。