Main Content

このページの内容は最新ではありません。最新版の英語を参照するには、ここをクリックします。

firpm

Parks-McClellan 最適 FIR フィルターの設計

説明

b = firpm(n,f,a) では、n 次の FIR フィルターの n+1 係数を含む行ベクトル b が返されます。結果として得られるフィルターの周波数と振幅の特性は、ベクトル fa によって与えられる特性と一致します。

b = firpm(n,f,a,w)w を使用して周波数ビンを重み付けします。

b = firpm(n,f,a,ftype)'ftype' によって指定されるフィルター タイプを使用します。

b = firpm(n,f,a,lgrid) は、整数 lgrid を使用して周波数グリッドの密度を制御します。

[b,err] = firpm(___) では、err の最大リップル高が返されます。前の入力構文のいずれでもこれを使用できます。

[b,err,res] = firpm(___) では、周波数応答の特性が構造体 res として返されます。

b = firpm(n,f,fresp,w) では、FIR フィルターが返されます。この周波数-振幅特性は、関数ハンドル fresp が返す応答を最もよく近似します。

b = firpm(n,f,fresp,w,ftype) では、反対称 (奇数) フィルターが設計されます。ここで、ftype は、フィルターを微分器またはヒルベルト変換器として指定します。fresp を指定しない場合、既定の対称性を決定するために ftype の呼び出しが行われます。

すべて折りたたむ

Parks-McClellan アルゴリズムを使用して 17 次の FIR バンドパス フィルターを設計します。正規化された阻止帯域周波数 0.3π および 0.7π ラジアン/サンプルと正規化された通過帯域周波数 0.4π および 0.6π ラジアン/サンプルを指定します。目的の周波数応答と実際の周波数応答をプロットします。

f = [0 0.3 0.4 0.6 0.7 1];
a = [0 0 1 1 0 0];
b = firpm(17,f,a);

[h,w] = freqz(b,1,512);
plot(f,a,w/pi,abs(h))
legend('Ideal','firpm Design')
xlabel 'Radian Frequency (\omega/\pi)', ylabel 'Magnitude'

Figure contains an axes object. The axes object with xlabel Radian Frequency ( omega / pi ), ylabel Magnitude contains 2 objects of type line. These objects represent Ideal, firpm Design.

1,500 Hz の通過帯域カットオフ周波数および 2,000 Hz の阻止帯域カットオフ周波数をもつローパス フィルターを設計します。サンプリング周波数を 8,000 Hz に指定します。最大阻止帯域振幅を 0.01、最大通過帯域誤差 (リップル) を 0.001 に指定する必要があります。firpmord を使用して必要なフィルター次数、正規化周波数帯域エッジ、周波数帯域振幅、および重みを取得します。

[n,fo,ao,w] = firpmord([1500 2000],[1 0],[0.001 0.01],8000);
b = firpm(n,fo,ao,w);
fvtool(b,1)

Figure Figure 1: Magnitude Response (dB) contains an axes object. The axes object with title Magnitude Response (dB), xlabel Normalized Frequency ( times pi blank rad/sample), ylabel Magnitude (dB) contains an object of type line.

Parks-McClellan アルゴリズムを使用して、1 kHz でサンプリングされた信号で使用される 50 次等リップルの FIR バンドパス フィルターを作成します。

N = 50;
Fs = 1e3;

通過帯域の周波数の範囲を 200 ~ 300 Hz で指定し、通過帯域のいずれかの側の遷移領域の幅を 50 Hz に指定します。

Fstop1 = 150;
Fpass1 = 200;
Fpass2 = 300;
Fstop2 = 350;

最適化適合が、低周波数阻止帯域に 3 の重み、通過帯域に 1 の重み、および高周波数阻止帯域に 100 の重みを付けるようにフィルターを設計します。フィルターの振幅応答を表示します。

Wstop1 = 3;
Wpass = 1;
Wstop2 = 100;

b = firpm(N,[0 Fstop1 Fpass1 Fpass2 Fstop2 Fs/2]/(Fs/2), ...
    [0 0 1 1 0 0],[Wstop1 Wpass Wstop2]);

fvtool(b,1)

Figure Figure 1: Magnitude Response (dB) contains an axes object. The axes object with title Magnitude Response (dB), xlabel Normalized Frequency ( times pi blank rad/sample), ylabel Magnitude (dB) contains an object of type line.

入力引数

すべて折りたたむ

フィルター次数。正の実数スカラーとして指定します。

正規化周波数点。実数値のベクトルとして指定します。引数は、[0, 1] の範囲でなければなりません。ここで、1 はナイキスト周波数に対応します。ベクトルの要素数は、常に 2 の倍数になります。周波数は増加する順でなければなりません。

f で指定した点での目的の振幅。ベクトルとして指定します。fa は同じ長さでなければなりません。長さは偶数でなければなりません。

  • k が偶数の場合 (f(k),f(k+1)) 点間の周波数における希望する振幅は、点 (f(k),a(k)) と点 (f(k+1),a(k+1)) とを結ぶ線分です。

  • k が奇数の場合、点のペア (f(k),f(k+1)) 間の周波数において希望する振幅は指定されません。このような点の間の領域は、遷移領域または特定のアプリケーションにとって重要ではない領域です。

各周波数帯域の適合の調整に使用される重み付け。実数値ベクトルとして指定します。w の長さは、f および a の長さの半分になるため、帯域ごとに厳密に 1 つの重みが存在します。

奇数の対称な線形位相フィルター (III 型および IV 型) のフィルター タイプ。'hilbert' または 'differentiator' として指定します。

  • 'hilbert'b の出力係数は、関係式 b(k) = -b(n + 2 - k), k = 1, ..., n + 1 に従います。このクラスのフィルターにはヒルベルト変換器が含まれ、これは帯域全体で目的の振幅 1 をもちます。

    次に例を示します。

    h = firpm(30,[0.1 0.9],[1 1],'hilbert');
    

    上の式では、長さ 31 の近似 FIR ヒルベルト変換器が設計されます。

  • 'differentiator' — ゼロ以外の振幅帯域に対しては、低周波数での誤差が高周波数での誤差よりも大幅に小さくなるように、フィルターによって誤差に 1/f 倍の重みが付けられます。周波数に比例した振幅特性をもつ FIR 微分器に対しては、これらのフィルターでは最大相対誤差 (希望する振幅に対する最大誤差率) が最小化されます。

周波数グリッドの密度を制御します。このグリッドにはおよそ (lgrid*n)/(2*bw) の周波数点があります。ここで、bwf でカバーされる周波数帯域区間 [0,1] 全体の分数です。lgrid を増やすと、多くの場合等リップルフィルターにほぼ正確に一致するフィルターとなりますが、計算に時間がかかります。既定値 16 は、lgrid に指定できる最小値です。

周波数応答。関数ハンドルとして指定します。関数は、以下の構文によって firpm の内部から呼び出されます。

[dh,dw] = fresp(n,f,gf,w)

引数は firpm の引数と同様です。

  • n はフィルター次数です。

  • f は、0 ~ 1 (1 はナイキスト周波数) で単調増加する、正規化周波数帯域エッジのベクトルです。

  • gf は、firpm により指定された各周波数帯域に対し線形に内挿されたグリッド点のベクトルです。gf は応答関数が評価される周波数グリッドを決定し、opt 構造の fgrid フィールドにcfirpmによって返されたデータと同じデータを含みます。

  • w は、最適化中に使用される帯域ごとの正で実数の重みのベクトルです。w は、firpm を呼び出す際にはオプションです。指定しない場合は、fresp に渡される前に重み 1 に設定されます。

  • dhdw は、それぞれ望ましい複素周波数応答と帯域重みのベクトルで、グリッド gf の各周波数で評価されます。

出力引数

すべて折りたたむ

フィルター係数。長さ n + 1 の行ベクトルとして返されます。係数は増加する順です。

最大リップル高。スカラーとして返されます。

周波数応答の特性。構造体として返されます。構造体 res には次のフィールドが含まれています。

res.fgrid

フィルター設計の最適化に使用される周波数グリッド ベクトル

res.des

res.fgrid の各点に対する望ましい周波数応答

res.wt

opt.fgrid の各点の重み

res.H

res.fgrid の各点に対する実際の周波数応答

res.error

res.fgrid (res.des-res.H) の各点の誤差

res.iextr

極値周波数に対する res.fgrid へのインデックスのベクトル

res.fextr

極値周波数のベクトル

ヒント

フィルター設計が収束しない場合、フィルター設計が正しくない可能性があります。周波数応答を調べて設計を検証してください。

フィルター設計が収束せず、生成されたフィルター設計が正しくない場合は、以下のうち 1 つまたは複数を試してください。

  • フィルター次数を大きくします。

  • 阻止帯域の減衰を小さくするか遷移領域を拡大する (あるいはその両方を行う) ことで、フィルターの設計を緩和します。

アルゴリズム

firpm では、Parks-McClellan アルゴリズム [2] を使用して、線形位相 FIR フィルターが設計されます。Parks-McClellan アルゴリズムは、Remez 交換アルゴリズムとチェビシェフ近似理論を使用し、希望する周波数応答と実際の周波数応答が最も近づくようにフィルターを設計します。フィルターは、希望周波数応答と実際の周波数応答間の最大誤差が最小化されているという点で、最適化されているといえます。この方法で設計されたフィルターは、周波数応答で等リップル動作を示すため、等リップル フィルターと呼ばれる場合もあります。firpm は、この等リップル特性のため、インパルス応答の最初と最後に不連続性を示します。

これらは、I 型 (次数 n が奇数) と II 型 (n が偶数) の線形位相フィルターになります。ベクトル fa では、フィルターの周波数-振幅特性が指定されます。

  • f は、0 ~ 1 の範囲で指定される周波数点のベクトルで、1 はナイキスト周波数に対応します。周波数は増加する順でなければなりません。

  • a は、f で指定した点で目的の振幅を含むベクトルです。

    k が奇数の場合、(f(k), f(k+1)) 両点間の周波数における目的の振幅関数は、点 (f(k), a(k)) と点 (f(k+1), a(k+1)) を結ぶ線分です。

    k が偶数の場合、(f(k), f(k+1)) 両点間の周波数における目的の振幅関数は指定されません。これらは遷移領域、つまり "don't care" 領域です。

  • fa は、同じ長さです。この長さは、偶数でなければなりません。

目的の振幅応答を定義する際のベクトル fa の関係は、以下のようになります。

firpm では、偶数の対称性と、ナイキスト周波数で非ゼロの通過帯域をもつ構成に対し、常に偶数のフィルター次数が使用されます。これは、偶数のフィルター次数が偶数の対称性と奇数の次数を示すインパルス応答に対し、ナイキスト周波数での周波数応答が必ず 0 になるからです。n に奇数の値を指定した場合、firpm によりこの値に 1 が加算されます。

firpm では、I 型、II 型、III 型、および IV 型の線形位相フィルターが設計されます。I 型および II 型は、それぞれ偶数 n および奇数 n の既定フィルターです。III 型 (n が偶数) および IV 型 (n が奇数) は、ftype 引数を使用して、それぞれ 'hilbert' または 'differentiator' で指定されます。異なるタイプのフィルターは、周波数応答に対して異なる対称性および制約をもちます (詳細については、[3]を参照してください。)

線形位相フィルター タイプフィルター次数係数の対称性f = 0 での応答 H(f)f = 1 (ナイキスト) での応答 H(f)

I 型

偶数

偶数:

b(k)=b(n+2k),k=1,...,n+1

制約なし

制約なし

II 型

奇数

偶数:

b(k)=b(n+2k),k=1,...,n+1

制約なし

H(1) = 0

firpm では、ナイキスト周波数で非ゼロの通過帯域をもつ II 型のフィルターを作成する場合、フィルターの次数を 1 増やします。

III 型

偶数

奇数:

b(k)=b(n+2k),k=1,...,n+1

H(0) = 0

H(1) = 0

IV 型奇数

奇数:

b(k)=b(n+2k),k=1,...,n+1

H(0) = 0

制約なし

firpm を使用して、希望する周波数応答を定義する関数を作成することもできます。firpm に対してあらかじめ定義された周波数応答関数のハンドルは @firpmfrf で、線形位相 FIR フィルターを設計するものです。

メモ

b = firpm(n,f,a,w) は、b = firpm(n,f,{@firpmfrf,a},w) と等価です。ここで、@firpmfrffirpm に対してあらかじめ定義された周波数応答関数ハンドルです。必要に応じて、自分独自の応答関数を作成できます。詳細については、help private/firpmfrf を使用して関数ハンドルの作成を参照してください。

参照

[1] Digital Signal Processing Committee of the IEEE Acoustics, Speech, and Signal Processing Society, eds. Selected Papers in Digital Signal Processing. Vol. II. New York: IEEE Press, 1976.

[2] Digital Signal Processing Committee of the IEEE Acoustics, Speech, and Signal Processing Society, eds. Programs for Digital Signal Processing. New York: IEEE Press, 1979, algorithm 5.1.

[3] Oppenheim, Alan V., Ronald W. Schafer, and John R. Buck. Discrete-Time Signal Processing. Upper Saddle River, NJ: Prentice Hall, 1999, p. 486.

[4] Parks, Thomas W., and C. Sidney Burrus. Digital Filter Design. New York: John Wiley & Sons, 1987, p. 83.

[5] Rabiner, Lawrence R., James H. McClellan, and Thomas W. Parks. "FIR Digital Filter Design Techniques Using Weighted Chebyshev Approximation." Proceedings of the IEEE®. Vol. 63, Number 4, 1975, pp. 595–610.

拡張機能

バージョン履歴

R2006a より前に導入