В этом примере показано, как проектировать пиковые и режущие фильтры. Фильтры, которые имеют пик или пробку на определенной частоте, полезны для сохранения или исключения конкретной частотной составляющей сигнала. Конструктивными параметрами фильтра являются частота, при которой требуется пик или пробка, и либо ширина полосы пропускания 3-dB, либо Q-фактор фильтра. Кроме того, учитывая эти характеристики, путем увеличения порядка фильтрации можно получить конструкции, которые более близко приближаются к идеальному фильтру.
Предположим, что необходимо устранить интерференцию 60 Гц в сигнале, дискретизированном на частоте 3000 Гц. Для этой цели можно использовать фильтр-насечку.
F0 = 60; % interference is at 60 Hz Fs = 3000; % sampling frequency is 3000 Hz notchspec = fdesign.notch('N,F0,Q',2,F0,10,Fs); notchfilt = design(notchspec,'SystemObject',true); fvtool(notchfilt,'Color','white');
Коэффициент качества или Q-фактор фильтра является мерой того, насколько хорошо желаемая частота изолирована от других частот. Для фиксированного порядка фильтрации более высокий Q-коэффициент достигается путем проталкивания полюсов ближе к нулям.
notchspec.Q = 100; notchfilt1 = design(notchspec,'SystemObject',true); fvt= fvtool(notchfilt, notchfilt1, 'Color','white'); legend(fvt,'Q = 10','Q = 100');
Эквивалентный способ определения коэффициента качества для определения полосы пропускания 3-dB, BW. Они связаны по Q = F0/BW. Определение полосы пропускания может быть более удобным способом достижения точно желаемой формы для проектируемого фильтра.
notchspec = fdesign.notch('N,F0,BW',2,60,5,3000); notchfilt2 = design(notchspec,'SystemObject',true); fvt= fvtool(notchfilt, notchfilt1, notchfilt2, 'Color','white'); legend(fvt,'Q = 10','Q = 100','BW = 5 Hz');
Поскольку до сих пор можно только толкать полюса и оставаться стабильными, чтобы улучшить приближение фильтра к кирпичной стене, необходимо увеличить порядок фильтра.
notchspec = fdesign.notch('N,F0,Q',2,.4,100); notchfilt = design(notchspec,'SystemObject',true); notchspec.FilterOrder = 6; notchfilt1 = design(notchspec,'SystemObject',true); fvt= fvtool(notchfilt, notchfilt1, 'Color','white'); legend(fvt,'2nd-Order Filter','6th-Order Filter');
Для данного порядка мы можем получить более резкие переходы, допуская пульсации полосы пропускания и/или полосы останова.
N = 8; F0 = 0.4; BW = 0.1; notchspec = fdesign.notch('N,F0,BW',N,F0,BW); notchfilt = design(notchspec,'SystemObject',true); notchspec1 = fdesign.notch('N,F0,BW,Ap,Ast',N,F0,BW,0.5,60); notchfilt1 = design(notchspec1,'SystemObject',true); fvt= fvtool(notchfilt, notchfilt1, 'Color','white'); legend(fvt,'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; peakspec = fdesign.peak('N,F0,BW',N,F0,BW); peakfilt = design(peakspec,'SystemObject',true); peakspec1 = fdesign.peak('N,F0,BW,Ast',N,F0,BW,80); peakfilt1 = design(peakspec1,'SystemObject',true); fvt= fvtool(peakfilt, peakfilt1, 'Color','white'); legend(fvt,'Maximally Flat 6th-Order Filter',... '6th-Order Filter With 80 dB Stopband Attenuation','Location','North');
Использование изменяющихся во времени фильтров требует изменения коэффициентов фильтра во время выполнения моделирования. Для дополнения процесса автоматического проектирования фильтров на основе объектов fdesign, DSP System Toolbox предоставляет другие возможности, включая функции для непосредственного вычисления коэффициентов фильтра, например iirnotch
Полезной отправной точкой является статический фильтр, вызываемый из динамического (потокового) моделирования. В этом случае непосредственно создается фильтр надрезов 2-го порядка, и его коэффициенты вычисляются с помощью iirnotch. Расчетными параметрами являются центральная частота 1 кГц и полоса пропускания 50 Гц при -3 дБ с частотой дискретизации 8 кГц.
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)
b = 1×3
0.8341 -1.1796 0.8341
a = 1×3
1.0000 -1.1796 0.6682
biquad = dsp.BiquadFilter('SOSMatrix', [b, a]); scope = dsp.SpectrumAnalyzer('PlotAsTwoSidedSpectrum', false, ... 'SampleRate', Fs, ... 'SpectralAverages', 16); samplesPerFrame = 256; nFrames = 4096; for k = 1:nFrames x = randn(samplesPerFrame, 1); y = biquad(x); scope(y); end
Коэффициенты изменяющихся во времени фильтров должны изменяться с течением времени из-за изменений во времени выполнения проектных параметров (например, центральной частоты для фильтра-пробки). Для каждого изменения расчетных параметров необходимо повторно вычислить векторы коэффициентов b и биквад. Для параметра SOSMatrix установлено новое значение. Это может быть дорогостоящим в вычислительном отношении. Для этого типа приложений dsp. OAllpassFilter может предложить более удобные конструкции фильтров. Преимущества включают - Внутренняя устойчивость - Коэффициенты, разделенные по проектным параметрам
Построить связанный решетчатый фильтр allpass, эквивалентный бикваду
[k1, k2] = tf2cl(b, a) %#okk1 =
[]
k2 = 2×1
-0.7071
0.6682
apnotch = dsp.CoupledAllpassFilter('Lattice', k1, k2); fvtool(apnotch, 'Fs', Fs)
Одно из преимуществ этой структуры, основанной на allpass, состоит в том, что новые коэффициенты разъединяются относительно конструктивных параметров F0 и BW. Так, например, изменение F0 на 3 кГц даст
F0 = 3e3/(Fs/2); [b, a] = iirnotch(F0, BW)
b = 1×3
0.8341 1.1796 0.8341
a = 1×3
1.0000 1.1796 0.6682
[k1, k2] = tf2cl(b, a)%#okk1 =
[]
k2 = 2×1
0.7071
0.6682
ПРИМЕЧАНИЕ: в то время как а и b оба изменились, в форме allpass решетки изменение конструкции затронуло только k2 (1).
Теперь измените полосу пропускания на 1 кГц
BW = 1e3/(Fs/2); [b, a] = iirnotch(F0, BW)
b = 1×3
0.7071 1.0000 0.7071
a = 1×3
1.0000 1.0000 0.4142
[k1, k2] = tf2cl(b, a)%#okk1 =
[]
k2 = 2×1
0.7071
0.4142
Изменение конструкции теперь затронуло только k2 (2).
Разделение коэффициентов имеет многочисленные преимущества в системах реального времени, включая более экономичное обновление коэффициентов и более предсказуемое переходное поведение при изменении коэффициентов.
Ниже приведен принцип изменения проектных параметров при динамическом моделировании, включая живую визуализацию влияния на оценку передаточной функции фильтра. В реальных приложениях обычно можно идентифицировать отдельные выражения, которые связывают каждый параметр конструкции с соответствующим коэффициентом allpass решетки. Использование этих функций вместо функций для всего фильтра, таких как 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 apnotch = dsp.CoupledAllpassFilter('Lattice'); scope = dsp.SpectrumAnalyzer('PlotAsTwoSidedSpectrum', false, ... 'SampleRate', Fs, ... 'SpectralAverages', 4, ... 'RBWSource', 'Auto'); paramtbl = ParameterTimeTable('Time', paramsChangeTimes, ... 'Values', myChangingParams, ... 'SampleRate', Fs/nSamplesPerFrame); % Actual simulation loop for frameIdx = 1:nFrames % Get current F0 and BW [params, update] = paramtbl(); 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 apnotch.LatticeCoefficients1{1} = k1; apnotch.LatticeCoefficients2{1} = k2; end % Generate vector of white noise samples x = randn(nSamplesPerFrame, 1); % Filter noise y = apnotch(x); % Visualize spectrum scope(y); end
