Этот пример показов, как разработать эффективные конечные импульсные характеристики фильтры с очень узкими полосами переходов с помощью многоступенчатых методов. Методы могут быть распространены на проект многоскоростных фильтров. Смотрите Многоступенчатый Проект Дециматоров/Интерполяторов для примера этого.
Одним из недостатков использования конечных импульсных характеристик является то, что порядок фильтра имеет тенденцию расти обратно пропорционально переходной полосе фильтра. Рассмотрите следующие проектные спецификации (где рябь указана в линейных модулях):
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 достигает эффективного проекта для вышеуказанных спецификаций в том смысле, что он уменьшает общее количество требуемых умножителей. Для этого задача проекта разбита на два этапа, фильтр, который усиливается для достижения строгих спецификаций, не используя много умножителей, и фильтр, который удаляет изображения, созданные при увеличении дискретизации предыдущего фильтра.
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')
В предыдущем примере мы автоматически определили коэффициент увеличения дискретизации, используемый таким образом, чтобы общее количество умножителей было сведено к минимуму. Выяснилось, что для заданных спецификаций оптимальный коэффициент увеличения дискретизации составил 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')
Заметьте, что затухание в полосе задерживания для многоступенчатой конструкции примерно в два раза больше, чем для других проектов. Это связано с тем, что дециматоры должны ослаблять компоненты вне полосы частот на требуемые 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');
Теперь мы показываем на примере, что 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