実践に即したデジタル フィルター設計の紹介

FIR フィルターの設計

ローパス フィルターの仕様

理想的なローパス フィルターは、信号のすべての周波数成分のうち、指定されたカットオフ周波数 $${\omega }_c$$ より小さい成分をすべて残し、$${\omega }_c$$ より大きい成分をすべて棄却するフィルターです。理想的なローパス フィルターの実装に必要なインパルス応答は無限に長くなるため、理想的な FIR ローパス フィルターを設計することはできません。理想的なインパルス応答での有限長の近似を行うと、フィルターの通過帯域 ($$\omega <{\omega }_c$$) および阻止帯域 ($$\omega >{\omega }_c$$) の両方でリップルが発生するほか、通過帯域と阻止帯域の間の遷移幅が非ゼロになります。

通過帯域リップルおよび阻止帯域リップルと遷移幅の発生は望ましくありませんが、有限のインパルス応答で近似を行う場合には、理想的なローパスの応答からの偏差として避けることはできません。この偏差は以下の図で示されます。

フィルター設計での設計仕様の位置付けは、各仕様を下の図に示す三角形の頂点の 1 つと考えるとわかりやすくなります。

三角形は、設計仕様を選択する際に使用可能な自由度を示すために使用されています。角度の合計は一定であるため、選択できるのは 2 つの仕様までです。3 番目の仕様は該当する設計アルゴリズムによって決定されます。また、三角形での角度と同様に、仕様の 1 つを大きくするか小さくすると、他の仕様の一方または両方に影響します。

FIR フィルターは本質的に安定性が高く線形位相をもつように設計できるため、非常に優れています。ただし、このフィルターは過渡応答が長く、適用する場合によっては非常に多くの計算を必要とすることがあります。

最小次数の FIR 設計

最小次数設計を行うには、通過帯域周波数、阻止帯域周波数、通過帯域リップルおよび阻止帯域の減衰量を指定します。これにより、設計アルゴリズムが、仕様を満たす最小のフィルター長を選択します。

通過帯域周波数 0.37*pi ラジアン/サンプル、阻止帯域周波数 0.43*pi ラジアン/サンプル (よって、遷移幅は 0.06*pi ラジアン/サンプル)、通過帯域リップル 1 dB、阻止帯域の減衰量 30 dB をもつ最小次数ローパス FIR フィルターを設計します。

Fpass = 0.37;
Fstop = 0.43;
Ap = 1;
Ast = 30;
d = designfilt("lowpassfir", ...
    PassbandFrequency=Fpass, ...
    StopbandFrequency=Fstop, ...
    PassbandRipple=Ap, ...
    StopbandAttenuation=Ast);
fa = filterAnalyzer(d,FilterNames="Equiripple_Design");

設計された次数は関数 filtord を使用してクエリできます。

N = filtord(d)

N = 39

フィルター設計に使用するパラメーターの詳細については、関数 info を使用して取得できます。

info(d)

ans = 25x32 char array
    'FIR Digital Filter (real)       '
    '-------------------------       '
    'Filter Length  : 40             '
    'Stable         : Yes            '
    'Linear Phase   : Yes (Type 2)   '
    '                                '
    'Design Method Information       '
    'Design Algorithm : Equiripple   '
    '                                '
    'Design Options                  '
    'Density Factor : 16             '
    'Maximum Phase  : false          '
    'Minimum Order  : any            '
    'Minimum Phase  : false          '
    'Stopband Decay : 0              '
    'Stopband Shape : flat           '
    'Uniform Grid   : true           '
    '                                '
    'Design Specifications           '
    'Sample Rate     : 2 (normalized)'
    'Response        : Lowpass       '
    'Stopband Atten. : 30 dB         '
    'Passband Edge   : 0.37          '
    'Stopband Edge   : 0.43          '
    'Passband Ripple : 1 dB          '

しかし、最小次数設計は、カイザー ウィンドウを使用しても可能なことに留意してください。カイザー ウィンドウ法では同じ仕様に対してフィルター次数が大きくなりますが、そのアルゴリズムは計算コストが低く、設計仕様が厳しい場合でも収束問題が発生する可能性が小さくなります。収束問題が発生する可能性があるのは、フィルターの適用において、遷移幅を非常に狭くする必要がある場合か、阻止帯域の減衰量が非常に大きい場合です。

カイザー ウィンドウを使用して上記と同じ仕様のフィルターを設計し、応答を等リップル フィルターと比較します。

dk = designfilt("lowpassfir", ...
    PassbandFrequency=Fpass, ...
    StopbandFrequency=Fstop, ...
    PassbandRipple=Ap, ...
    StopbandAttenuation=Ast, ...
    DesignMethod="kaiserwin");
addFilters(fa,dk,FilterNames="Kaiser_Window_Design");

N = filtord(dk)

N = 52

Hz での周波数パラメーターの指定

フィルターが動作するサンプル レートが既知の場合は、サンプル レートと周波数を Hz で指定できます。2 kHz のサンプル レートに対する最小次数等リップル フィルターを再設計します。

Fpass = 370;
Fstop = 430;
Ap = 1;
Ast = 30;
Fs = 2000;
d = designfilt("lowpassfir", ...
    PassbandFrequency=Fpass, ...
    StopbandFrequency=Fstop, ...
    PassbandRipple=Ap, ...
    StopbandAttenuation=Ast, ...
    SampleRate=Fs);
filterAnalyzer(d)

固定次数、固定遷移幅

固定次数設計は計算負荷の影響を受けやすいフィルター適用や、フィルター係数の数に制限のあるフィルター適用を行う場合に役立ちます。通過帯域リップル/阻止帯域の減衰量の制御を犠牲にして遷移幅を固定するという選択肢もあります。

通過帯域周波数が 370 Hz、阻止帯域周波数が 430 Hz、サンプル レートが 2 kHz の 30 次ローパス FIR フィルターについて考えます。この仕様セットには、2 つの設計法が使用可能です。つまり、等リップルと最小二乗です。それぞれの方法でフィルターを設計し、結果を比較しましょう。

N = 30;
Fpass = 370;
Fstop = 430;
Fs = 2000;
% Design method defaults to "equiripple" when omitted
deq = designfilt("lowpassfir", ...
    FilterOrder=N, ...
    PassbandFrequency=Fpass, ...
    StopbandFrequency=Fstop, ...
    SampleRate=Fs);
dls = designfilt("lowpassfir", ...
    FilterOrder=N, ...
    PassbandFrequency=Fpass, ...
    StopbandFrequency=Fstop, ...
    DesignMethod="ls", ...
    SampleRate=Fs);
filterAnalyzer(deq,dls,FilterNames=["Equiripple_Design_2" "Least_Squares_Design_2"])

等リップル フィルターは、指定された最小の阻止帯域の減衰量や指定された最大通過帯域リップルのあるフィルターを設計する場合など、特定の許容誤差を満たす必要があるフィルター適用に最適です。一方、この設計は、通過帯域や阻止帯域で (理想的なフィルターと実際のフィルターとの間の) 誤差のエネルギーを最小化する必要がある場合、望ましくない可能性があります。

前述の例では、設計されたフィルターには通過帯域と阻止帯域で同じリップルがありました。フィルター次数を固定しながら、重みを使用してどちらかの帯域でのリップルを減らすことができます。たとえば、阻止帯域リップルを通過帯域リップルの 1/10 にするには、阻止帯域に通過帯域の 10 倍の重みを指定しなければなりません。この事実を使用して等リップル フィルターを再設計します。

deqw = designfilt("lowpassfir", ...
    FilterOrder=N, ...
    PassbandFrequency=Fpass, ...
    StopbandFrequency=Fstop, ...
    PassbandWeight=1, ...
    StopbandWeight=10, ...
    SampleRate=Fs);
filterAnalyzer(deq,deqw,FilterNames=["Equiripple_Design_3" "Equiripple_Design_With_Weighted_Stopband"])

固定次数、固定カットオフ周波数

ウィンドウ設計法を使用して、固定フィルター次数および固定カットオフ周波数をもつフィルターを設計できます。

たとえば、カットオフ周波数が 60 Hz、サンプル レートが 1 kHz の 100 次ローパス FIR フィルターについて考えます。ハミング ウィンドウを使用して作成された設計と、チェビシェフ ウィンドウおよび 90 dB のサイド ローブ減衰を使用して作成された設計を比較します。

dhamming = designfilt("lowpassfir", ...
    FilterOrder=100, ...
    CutoffFrequency=60, ...
    SampleRate=1000, ...
    Window="hamming");
dchebwin = designfilt("lowpassfir", ...
    FilterOrder=100, ...
    CutoffFrequency=60, ...
    SampleRate=1000, ...
    Window={"chebwin",90});
filterAnalyzer(dhamming,dchebwin,FilterNames=["Hamming_Window" "Chebyshev_Window"])

固定次数のあるフィルターを指定するには、他の方法もあります。つまり、固定カットオフ周波数、通過帯域リップルおよび阻止帯域の減衰量によるもの、固定遷移幅によるもの、そして固定電力半値 (3 dB) 周波数によるものです。

IIR フィルターの設計

FIR フィルターの欠点の 1 つは、特定の設計仕様を満たすために大きなフィルター次数を必要とすることです。リップルが一定に維持されると、フィルター次数は遷移幅に逆比例して増大します。フィードバックを使用することで、はるかに小さなフィルター次数で一連の設計仕様を満たすことができます。これが IIR フィルター設計の背景にある概念です。"無限インパルス応答" (IIR) という用語は、そのフィルターにインパルスが適用された場合に、出力がの減衰がゼロにならないという事実に由来しています。

IIR フィルターを使用するもう 1 つの重要な理由は、FIR フィルターに比べて群遅延が小さいため、過渡特性応答が短くて済むためです。

バタワース フィルター

バタワース フィルターは最大フラット IIR フィルターです。通過帯域と阻止帯域が平坦であると、遷移帯域が非常に広くなります。狭い遷移幅のフィルターを得るには、大きな次数が必要になります。

通過帯域周波数が 100 Hz、阻止帯域周波数が 300 Hz、最大通過帯域リップルが 1 dB、阻止帯域の減衰量が 60 dB である最小次数バタワース フィルターを設計します。サンプル レートは 2 kHz です。

Fp = 100;
Fst = 300;
Ap = 1;
Ast = 60;
Fs = 2e3;
dbutter = designfilt("lowpassiir", ...
    PassbandFrequency=Fp, ...
    StopbandFrequency=Fst, ...
    PassbandRipple=Ap, ...
    StopbandAttenuation=Ast, ...
    SampleRate=Fs, ...
    DesignMethod="butter");

チェビシェフ I 型フィルター

チェビシェフ I 型フィルターは、通過帯域リップルを許容することで、同じ次数のバタワース フィルターより遷移幅が小さくなります。

前述のバタワース フィルターと同じ仕様のチェビシェフ I 型フィルターを設計します。

dcheby1 = designfilt("lowpassiir", ...
    PassbandFrequency=Fp, ...
    StopbandFrequency=Fst, ...
    PassbandRipple=Ap, ...
    StopbandAttenuation=Ast, ...
    SampleRate=Fs, ...
    DesignMethod="cheby1");

チェビシェフ II 型フィルター

極端に大きな減衰量は通常必要ないため、いずれかの阻止帯域リップルを許容することで比較的小さな次数で必要な遷移幅を得ることができる可能性があります。

前の例と同じ仕様の最小次数チェビシェフ II 型フィルターを設計します。

dcheby2 = designfilt("lowpassiir", ...
    PassbandFrequency=Fp, ...
    StopbandFrequency=Fst, ...
    PassbandRipple=Ap, ...
    StopbandAttenuation=Ast, ...
    DesignMethod="cheby2", ...
    SampleRate=Fs);

楕円フィルター

楕円フィルターは、通過帯域と阻止帯域の両方でリップルを許容することで、チェビシェフ フィルターとバタワース フィルターを汎用的にします。リップルを小さくすると、楕円フィルターは、チェビシェフ フィルターまたはバタワース フィルターの位相応答および振幅を任意の近傍で近似することができます。

dellip = designfilt("lowpassiir", ...
    PassbandFrequency=Fp, ...
    StopbandFrequency=Fst, ...
    PassbandRipple=Ap, ...
    StopbandAttenuation=Ast, ...
    DesignMethod="ellip", ...
    SampleRate=Fs);

4 つの IIR フィルターの応答と次数を比較します。

FilterOrders = [filtord(dbutter) filtord(dcheby1) filtord(dcheby2) filtord(dellip)]

FilterOrders = 1×4

     7     5     5     4

fa = filterAnalyzer(dbutter,dcheby1,dcheby2,dellip, ...
    FilterNames=["Butterworth","Chebyshev_Type_I" "Chebyshev_Type_II" "Elliptic"]);
zoom(fa,"xy",[0 1e3 -80 2])

通過帯域を拡大し、リップルの差を確認します。

zoom(fa,"xy",[0 150 -3 2])

通過帯域または阻止帯域仕様への完全な準拠

最小次数設計では、理想的な次数を次の整数に丸めなければなりません。この追加の非整数次数により、アルゴリズムは実際に仕様を超えることになります。

既定の設定では、チェビシェフ I 型設計は通過帯域に一致します。バタワースとチェビシェフ II 型は阻止帯域に一致し、楕円設計は通過帯域と阻止帯域の両方に一致します (阻止帯域のエッジ周波数は超えます)。

dellip1 = designfilt("lowpassiir", ...
    PassbandFrequency=Fp, ...
    StopbandFrequency=Fst, ...
    PassbandRipple=Ap, ...
    StopbandAttenuation=Ast, ...
    DesignMethod="ellip", ...
    SampleRate=Fs, ...
    MatchExactly="passband");

dellip2 = designfilt("lowpassiir", ...
    PassbandFrequency=Fp, ...
    StopbandFrequency=Fst, ...
    PassbandRipple=Ap, ...
    StopbandAttenuation=Ast, ...
    DesignMethod="ellip", ...
    SampleRate=Fs, ...
    MatchExactly="stopband");
fa = filterAnalyzer(dellip,dellip1,dellip2, ...
    FilterNames=["Matched_Passband_And_Stopband" ...
    "Matched_Passband" "Matched_Stopband"]);

zoom(fa,"xy",[0 1e3 -80 2])

一致する通過帯域と一致する両方の設計のリップルは、通過帯域周波数値 100 Hz で正確に 1 dB です。

群遅延の比較

IIR フィルターの場合、リップルと遷移幅のトレードオフだけでなく位相歪みの度合いも考慮する必要があります。ナイキスト区間全体にわたって線形位相をもつことが不可能であることはわかっています。このため、位相応答がどの程度線形でないかを確認する必要があります。これを行うには、(理想的には一定の) 群遅延を見てそれがどの程度平坦であるかを確認します。

前に設計した 4 つの IIR フィルターの群遅延を比較します。

fa = filterAnalyzer(dbutter,dcheby1,dcheby2,dellip,Analysis="groupdelay", ...
    FilterNames=["Butterworth_1" "Chebyshev_Type_I_2" "Chebyshev_Type_II_2" "Elliptic_2"]);

まとめ

この例では、designfilt を使用して、さまざまな制約と設計法をもついろいろなローパス FIR および IIR フィルターを作成する方法を学習しました。designfilt は、ハイパス、バンドパス、バンドストップ、任意振幅、微分器およびヒルベルトの設計にも使用できます。すべての利用可能なオプションの詳細については、フィルター設計ギャラリーを参照してください。

参考情報

フィルター設計および解析についての詳細については、Signal Processing Toolbox® ソフトウェア ドキュメンテーションを参照してください。フィルター適用についての詳細については、実践に即したデジタル フィルター処理の紹介の例を参照してください。

Copyright 2024 The MathWorks, Inc.