ドキュメンテーション センター

  • 評価版
  • 製品アップデート

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

ピーク フィルターおよびノッチ フィルターの設計

この例では、ピーク フィルターおよびノッチ フィルターの設計法を示します。特定の周波数でピークに達するまたは刻むフィルターは、信号の特定の周波数コンポーネントを保持または削除する場合に便利です。フィルターの設計パラメーターは、ピークまたはノッチに指定された周波数、および 3-dB の帯域幅またはフィルターの Q ファクターのいずれかです。さらに、指定されたこれらの仕様では、フィルター次数を増加することで、理想フィルターに密接に近似する設計を得ることが可能です。

ノッチ フィルター

3000 Hz で抽出された信号の 60 Hz の干渉を排除する必要があるとします。このような目的にはノッチ フィルターを使用できます。

F0 = 60;   % interference is at 60 Hz
Fs = 3000; % sampling frequency is 3000 Hz
f = fdesign.notch('N,F0,Q',2,F0,10,Fs);
h = design(f,'SystemObject',true);
hfvt= fvtool(h,'Color','white');

フィルターのクオリティ ファクターまたは Q ファクターは、希望の周波数における他の周波数からの分離性能を示す尺度です。固定フィルター次数の場合、高次の Q ファクターは極を押してゼロに近づけることで実行されます。

f.Q = 100;
h1 = design(f,'SystemObject',true);
hfvt= fvtool(h, h1, 'Color','white');
legend(hfvt,'Q = 10','Q = 100');

クオリティ ファクターを指定する他の方法としては、3-dB の帯域幅 BW を指定する方法があります。これらは Q = F0/BW に関連しています。帯域幅を指定する方が、設計されたフィルターの最適な形状を実現するために便利な方法である場合があります。

f = fdesign.notch('N,F0,BW',2,60,5,3000);
h2 = design(f,'SystemObject',true);
hfvt= fvtool(h, h1, h2, 'Color','white');
legend(hfvt,'Q = 10','Q = 100','BW = 5 Hz');

これまで極を押して安定させることのみが可能だったので、フィルターの急峻な特性をもつ近似を改良するためには、フィルターの次数を増加させる必要があります。

f = fdesign.notch('N,F0,Q',2,.4,100);
h = design(f,'SystemObject',true);
f.FilterOrder = 6;
h1 = design(f,'SystemObject',true);
hfvt= fvtool(h, h1, 'Color','white');
legend(hfvt,'2nd-Order Filter','6th-Order Filter');

指定された次数では、通過帯域および/または阻止帯域リップルを可能にすることで鋭い遷移を得ることができます。

N = 8; F0 = 0.4; BW = 0.1;
f = fdesign.notch('N,F0,BW',N,F0,BW);
h = design(f,'SystemObject',true);
f1 = fdesign.notch('N,F0,BW,Ap,Ast',N,F0,BW,0.5,60);
h1 = design(f1,'SystemObject',true);
hfvt= fvtool(h, h1, 'Color','white');
legend(hfvt,'Maximally Flat 8th-Order Filter',...
    '8th-Order Filter With Passband/Stopband Ripples', ...
    'Location','NorthWest');
axis([0 1 -90 0.5]);

ピーク フィルター

ピーク フィルターは、信号から単一の周波数成分 (または周波数の小さな帯域) のみを保持する場合に使用されます。これまでに説明したすべての仕様およびトレードオフは、ピーク フィルターに同様に適用されます。例を以下に示します。

N = 6; F0 = 0.7; BW = 0.001;
f = fdesign.peak('N,F0,BW',N,F0,BW);
h = design(f,'SystemObject',true);
f1 = fdesign.peak('N,F0,BW,Ast',N,F0,BW,80);
h1 = design(f1,'SystemObject',true);
hfvt= fvtool(h, h1, 'Color','white');
legend(hfvt,'Maximally Flat 6th-Order Filter',...
    '6th-Order Filter With 80 dB Stopband Attenuation','Location','North');

時変ノッチ フィルターの実装

時変フィルターを使用するには、シミュレーションの実行中にフィルターの係数を変更しなければなりません。fdesign オブジェクトに基づく自動フィルター設計ワークフローを補完するために、DSP System Toolbox ではその他の機能が提供されており、iirnotch などのフィルター係数を直接計算する関数も含まれています。

動的な (ストリーミングされた) シミュレーション内から呼び出される静的フィルターから始めるとしましょう。この場合、2 次ノッチ フィルターが直接作成され、iirnotch で係数が計算されます。設計パラメーターは -3 dB において中心周波数が 1 kHz、帯域幅が 50 Hz、サンプリング周波数は 8 kHz です。

Fs = 8e3;           % 8 kHz sampling frequency
F0 = 1e3/(Fs/2);    % Notch at 2 kHz
BW = 500/(Fs/2);    % 500 Hz bandwidth
[b, a] = iirnotch(F0, BW)
hbq = dsp.BiquadFilter('SOSMatrix', [b, a]);

hspectr = dsp.SpectrumAnalyzer('PlotAsTwoSidedSpectrum', false, ...
    'SampleRate', Fs, ...
    'SpectralAverages', 16);
samplesPerFrame = 256;

nFrames = 4096;
for k = 1:nFrames
   x = randn(samplesPerFrame, 1);
   y = step(hbq, x);
   step(hspectr, y);
end
b =

    0.8341   -1.1796    0.8341


a =

    1.0000   -1.1796    0.6682

設計パラメーター (ノッチ フィルターの中心周波数など) は実行時に変化するため、時変フィルターの係数は時間経過とともに変化しなければなりません。設計パラメーターが変化するたびに、係数ベクトル b と a を再計算し、hbq.SOSMatrix を新しい値に設定しなければなりません。これには大量の計算が必要です。このようなタイプのアプリケーションでは、dsp.CoupledAllpassFilter の方が便利なフィルター構造を提供する場合があります。利点としては、本質的な安定性と、係数の設計パラメーターからの分離があげられます。

% Build a coupled allpass lattice filter equivalent to hbq
[k1, k2] = tf2cl(b, a)
hnotch = dsp.CoupledAllpassFilter('Lattice', k1, k2);

fvtool(hnotch, 'Fs', Fs)
k1 =

     []


k2 =

   -0.7071
    0.6682

このオールパスベースの構造の利点としては、新しい係数が設計パラメーター F0 および BW から分離されているという点があります。このため、たとえば F0 を 3 kHz にすると、出力は次のようになります。

F0 = 3e3/(Fs/2);
[b, a] = iirnotch(F0, BW)
[k1, k2] = tf2cl(b, a)

% NOTE: while a and b both changed, in the lattice allpass form the design
% change only affected k2(1)

% Now change the bandwidth to 1 kHz
BW = 1e3/(Fs/2);
[b, a] = iirnotch(F0, BW)
[k1, k2] = tf2cl(b, a)

% The design change now only affected k2(2)

% Coefficient decoupling has numerous advantages in real-time systems,
% including more economical coefficient update and more predictable
% transient behaviour when coefficients change
b =

    0.8341    1.1796    0.8341


a =

    1.0000    1.1796    0.6682


k1 =

     []


k2 =

    0.7071
    0.6682


b =

    0.7071    1.0000    0.7071


a =

    1.0000    1.0000    0.4142


k1 =

     []


k2 =

    0.7071
    0.4142

以下は上記の原理を適用して、動的なシミュレーション中に設計パラメーターを変更しています。また、フィルター伝達関数の推定に対する影響をライブで可視化しています。実際のアプリケーションでは、一般的に各設計パラメーターを対応するラティス オールパス係数にリンクする個々の式を識別します。主シミュレーション ループで iirnotch や tf2cl などのフィルター全体の関数の代わりにこれらの式を使用することで、効率性が向上されます。

% Notch filter parameters - how they vary over time
Fs = 8e3;       % 8 kHz sampling frequency
f0 = 1e3*[0.5, 1.5, 3, 2]/(Fs/2);   % in [0, 1] range
bw = 500/(Fs/2) * ones(1,4);        % in [0, 1] range

myChangingParams = struct('F0', num2cell(f0), 'BW', num2cell(bw));
paramsChangeTimes = [0, 7, 15, 20];       % in seconds

% Simulation time management
nSamplesPerFrame = 256;
tEnd = 30;
nSamples = ceil(tEnd * Fs);
nFrames = floor(nSamples / nSamplesPerFrame);

% Object creation
hnotch = dsp.CoupledAllpassFilter('Lattice');
hspectr = dsp.SpectrumAnalyzer('PlotAsTwoSidedSpectrum', false, ...
    'SampleRate', Fs, ...
    'SpectralAverages', 4, ...
    'RBWSource', 'Auto');
htable = dspdemo.ParameterTimeTable('Time', paramsChangeTimes, ...
    'Values', myChangingParams, ...
    'SampleRate', Fs/nSamplesPerFrame);

% Actual simulation loop
for frameIdx = 1:nFrames
    % Get current F0 and BW
    [params, update] = step(htable);

    if(update)
        % Recompute filter coefficients if parameters changed
        [b, a] = iirnotch(params.F0, params.BW);
        [k1, k2] = tf2cl(b, a);
        % Set filter coefficients to new values
        hnotch.LatticeCoefficients1{1} = k1;
        hnotch.LatticeCoefficients2{1} = k2;
    end

    % Generate vector of white noise samples
    x = randn(nSamplesPerFrame, 1);
    % Filter noise
    y = step(hnotch, x);
    % Visualize spectrum
    step(hspectr, y);
end

この情報は役に立ちましたか?