Эффективная узкополосная переходная конечная импульсная характеристика Создания фильтра

Этот пример показов, как разработать эффективные конечные импульсные характеристики фильтры с очень узкими полосами переходов с помощью многоступенчатых методов. Методы могут быть распространены на проект многоскоростных фильтров. Смотрите Многоступенчатый Проект Дециматоров/Интерполяторов для примера этого.

Проект Lowpass с узкой переходной полосой

Одним из недостатков использования конечных импульсных характеристик является то, что порядок фильтра имеет тенденцию расти обратно пропорционально переходной полосе фильтра. Рассмотрите следующие проектные спецификации (где рябь указана в линейных модулях):

Fpass = 2866.5;  % Passband edge
Fstop = 3087;    % Stopband edge
Apass = 0.0174;  % Passband ripple, 0.0174 dB peak to peak
Astop = 66.0206; % Stopband ripple, 66.0206 dB of minimum attenuation
Fs    = 44.1e3;

lowpassSpec = fdesign.lowpass(Fpass,Fstop,Apass,Astop,Fs);

Регулярная линейно-фазовая конструкция равновесия, которая соответствует спецификациям, может быть спроектирована с:

eqripFilter = design(lowpassSpec,'equiripple','SystemObject',true);
cost(eqripFilter)
ans = struct with fields:
                  NumCoefficients: 695
                        NumStates: 694
    MultiplicationsPerInputSample: 695
          AdditionsPerInputSample: 694

Требуемая длина фильтра оказывается 694 отвода.

Проект интерполированной конечной импульсной характеристики (IFIR)

Алгоритм проекта IFIR достигает эффективного проекта для вышеуказанных спецификаций в том смысле, что он уменьшает общее количество требуемых умножителей. Для этого задача проекта разбита на два этапа, фильтр, который усиливается для достижения строгих спецификаций, не используя много умножителей, и фильтр, который удаляет изображения, созданные при увеличении дискретизации предыдущего фильтра.

interpFilter = design(lowpassSpec,'ifir','SystemObject',true);

Видимо, мы сделали всё хуже. Вместо одного фильтра с 694 умножителями у нас теперь есть два фильтра с общим числом 804 умножителей. Однако пристальное изучение второго этапа показывает, что только около одного множителя в 5 является ненулевым. Фактическое общее количество умножителей уменьшено с 694 до 208.

cost(interpFilter)
ans = struct with fields:
                  NumCoefficients: 208
                        NumStates: 802
    MultiplicationsPerInputSample: 208
          AdditionsPerInputSample: 206

Давайте сравним ответы двух проектов:

fvt = fvtool(eqripFilter,interpFilter,'Color','White');
legend(fvt,'Equiripple design', 'IFIR design','Location','Best')

Figure Filter Visualization Tool - Magnitude Response (dB) contains an axes and other objects of type uitoolbar, uimenu. The axes with title Magnitude Response (dB) contains 3 objects of type line. These objects represent Equiripple design, IFIR design.

Ручное управление коэффициентом повышающей дискретизации

В предыдущем примере мы автоматически определили коэффициент увеличения дискретизации, используемый таким образом, чтобы общее количество умножителей было сведено к минимуму. Выяснилось, что для заданных спецификаций оптимальный коэффициент увеличения дискретизации составил 5. Однако, если мы исследуем опции проекта:

opts = designopts(lowpassSpec,'ifir')
opts = struct with fields:
      FilterStructure: 'dffir'
     UpsamplingFactor: 'auto'
    JointOptimization: 0
         SystemObject: 0

мы видим, что мы можем контролировать коэффициент увеличения дискретизации. Для примера, если бы мы хотели улучшить на 4, а не на 5:

opts.UpsamplingFactor = 4;
opts.SystemObject = true; 
upfilter = design(lowpassSpec,'ifir',opts);
cost(upfilter)
ans = struct with fields:
                  NumCoefficients: 217
                        NumStates: 767
    MultiplicationsPerInputSample: 217
          AdditionsPerInputSample: 215

Мы получили бы проект, который имеет в общей сложности 217 ненулевых умножителей.

Использование оптимизации соединений

Возможно конъюнктурное проектирование двух фильтров, используемых в IFIR. Делая это, мы можем сэкономить значительное количество умножителей за счет более длительного времени проектирования (из-за особенностей алгоритма проект также может не сходиться вообще в некоторых случаях):

opts.UpsamplingFactor = 'auto'; % Automatically determine the best factor
opts.JointOptimization = true;
ifirFilter = design(lowpassSpec,'ifir',opts);
cost(ifirFilter)
ans = struct with fields:
                  NumCoefficients: 172
                        NumStates: 726
    MultiplicationsPerInputSample: 172
          AdditionsPerInputSample: 170

Для этого проекта был найден лучший коэффициент повышающей дискретизации 6. Количество ненулевых умножителей теперь только 152

Использование мультирейтовых/многоступенчатых методов для достижения эффективных проектов

Для рассмотренных до настоящего времени проектов использовались односкоростные методы. Это означает, что количество умножений, требуемых на входную выборку (MPIS), равно количеству ненулевых умножителей. Для образца последний проект, которую мы показали, требует 152 MPIS. Проект одноступенчатого равновесия, который мы начали с необходимого 694 MPIS.

Используя мультирейт/многоступенчатые методы, которые комбинируют десятикратное уменьшение и интерполяцию, мы также можем получить эффективные проекты с низким количеством MPIS. Для дециматоров количество умножений, требуемых на входную выборку (в среднем), определяется количеством умножителей, разделенных на коэффициент десятикратного уменьшения.

multistageFilter = design(lowpassSpec,'multistage','SystemObject',true)
multistageFilter = 
  dsp.FilterCascade with properties:

    Stage1: [1x1 dsp.FIRDecimator]
    Stage2: [1x1 dsp.FIRDecimator]
    Stage3: [1x1 dsp.FIRDecimator]
    Stage4: [1x1 dsp.FIRInterpolator]
    Stage5: [1x1 dsp.FIRInterpolator]
    Stage6: [1x1 dsp.FIRInterpolator]

cost(multistageFilter)
ans = struct with fields:
                  NumCoefficients: 396
                        NumStates: 352
    MultiplicationsPerInputSample: 73
          AdditionsPerInputSample: 70.8333

Первый каскад имеет 21 умножитель, и коэффициент десятикратного уменьшения 3. Поэтому количество MPIS составляет 7. Второй этап имеет длину 45 и совокупный коэффициент десятикратного уменьшения 6 (то есть коэффициент десятикратного уменьшения этого этапа, умноженный на коэффициент десятикратного уменьшения первого этапа; это связано с тем, что входные выборки к этому этапу уже поступают со скоростью 1/3 скорости входных выборок к первому этапу). Среднее количество умножений на входную выборку (ссылка на вход общего мультирейт/многоступенчатого фильтра), таким образом, 45/6 = 7,5. Наконец, учитывая, что третий каскад имеет коэффициент десятикратного уменьшения 1, среднее количество умножений на вход для этого каскада составляет 130/6 = 21,667. Общее количество средних MPIS для трех дециматоров составляет 36,5.

Для интерполяторов оказывается, что фильтры идентичны дециматорам. При этом их вычислительная стоимость одинаковая. Поэтому общее количество MPIS для всего многоуровневого / многоступенчатого проекта составляет 73.

Теперь мы сравниваем ответы проекта equiripple и этого:

fvt = fvtool(eqripFilter,multistageFilter,'Color','White', ...
    'NormalizeMagnitudeto1','on');
legend(fvt,'Equiripple design', 'Multirate/multistage design', ...
    'Location','NorthEast')

Figure Filter Visualization Tool - Magnitude Response (dB) contains an axes and other objects of type uitoolbar, uimenu. The axes with title Magnitude Response (dB) contains 3 objects of type line. These objects represent Equiripple design, Multirate/multistage design.

Заметьте, что затухание в полосе задерживания для многоступенчатой конструкции примерно в два раза больше, чем для других проектов. Это связано с тем, что дециматоры должны ослаблять компоненты вне полосы частот на требуемые 66 дБ, порядок избежать сглаживания, которое будет нарушать требуемые спецификации. Точно так же интерполяторам необходимо ослабить изображения на 66 дБ в порядок, чтобы удовлетворить спецификациям этой проблемы.

Управление количеством этапов вручную

Полученный многоступенчатый / многоступенчатый проект состоял из 6 этапов. Количество каскадов определяется автоматически по умолчанию. Однако также можно вручную управлять количеством этапов, которые приводят. Для примера:

lp4stage=design(lowpassSpec,'multistage','NStages',4,'SystemObject',true);
cost(lp4stage)
ans = struct with fields:
                  NumCoefficients: 516
                        NumStates: 402
    MultiplicationsPerInputSample: 86
          AdditionsPerInputSample: 84.5000

Среднее количество MPIS для этого случая составляет 86

Групповая задержка

Мы можем вычислить групповую задержку для каждого проекта. Заметьте, что многоуровневая/многоступенчатая конструкция вводит наибольшую задержку (это цена для оплаты менее вычислительно дорогого проекта). Проект IFIR вводит большую задержку, чем одноступенчатая конструкция equiripple, но меньшую, чем многоступенчатая/многоступенчатая конструкция.

fvt = fvtool(eqripFilter,interpFilter,multistageFilter,...
    'Analysis','grpdelay','Color','White');
legend(fvt, 'Equiripple design','IFIR design',...
    'Multirate/multistage design');

Figure Filter Visualization Tool - Group delay contains an axes and other objects of type uitoolbar, uimenu. The axes with title Group delay contains 3 objects of type line. These objects represent Equiripple design, IFIR design, Multirate/multistage design.

Фильтрация сигнала

Теперь мы показываем на примере, что IFIR и многоступенчатая/мультирейтовая проект работают сопоставимо с одноступенчатым проектом equiripple, требуя гораздо меньших расчетов. В порядок выполнения фильтрации для каскадов необходимо вызвать generateFilteringCode (ifirFilter) и generateFilteringCode (multistageFilter). Это уже сделано здесь в порядок создания HelperIFIRFilter и HelperMultiFIRFilter.

scope = dsp.SpectrumAnalyzer('SpectralAverages',50,'SampleRate',Fs,...
    'PlotAsTwoSidedSpectrum',false,'YLimits', [-90 10],...
    'ShowLegend',true','ChannelNames',{'Equiripple design',...
    'IFIR design','Multirate/multistage design'});
tic,
while toc < 20
    % Run for 20 seconds
    x = randn(6000,1);
    
    % Filter using single stage FIR filter
    y1 = eqripFilter(x);
    
    % Filter using IFIR filter
    y2 = HelperIFIRFilter(x);
    
    % Filter multi-stage/multi-rate FIR filters
    y3 = HelperMultiFIRFilter(x);
    
    % Compare the output from both approaches
    scope([y1,y2,y3])
end

Figure Spectrum Analyzer contains an axes and other objects of type uiflowcontainer, uimenu, uitoolbar. The axes contains 3 objects of type line. These objects represent Equiripple design, IFIR design, Multirate/multistage design.