В этом примере показов, как разработать конечную импульсную характеристику полуполосы фильтров. Полуполосные фильтры широко используются в многочастотных приложениях обработки сигналов при интерполяции/децимировании в два раза. Полуполосные фильтры реализованы эффективно в полифазной форме, потому что примерно половина его коэффициентов равны нулю.
Полуполосы фильтры имеют две важные характеристики, ширина полосы пропускания и пульсация полосы пропускания должны быть одинаковыми, и частоты края полосы пропускания и края полосы остановки равноудалены от частоты полуполосы Fs/4 (Or pi/2 рад/отсчета в нормализованной частоте).
The firhalfband функция возвращает коэффициенты фильтра конечной импульсной характеристики полуполосы equiripple. В качестве простого примера рассмотрим полуполосный фильтр, чья работа с данными, дискретизированными с частотой 96 кГц и частотой полосы пропускания 22 кГц.
Fs = 96e3; Fp = 22e3; N = 100; num = firhalfband(N,Fp/(Fs/2)); fvt = fvtool(num,'Fs',Fs,'Color','white'); fvt.MagnitudeDisplay = 'Zero-phase';
Путем изменения масштаба отклика можно проверить, что полоса пропускания и полоса остановки от пика до пика ряби одинаковы. Также существует симметрия относительно Fs/4 (24 кГц) точки. Ширина полосы пропускания составляет до 22 кГц, как указано, и полоса остановки начинается с 26 кГц. Мы также можем проверить, что каждый другой коэффициент равен нулю, рассматривая импульсную характеристику. Это делает фильтр очень эффективным для реализации для интерполяции/децимации в 2 раза.
fvt.Analysis = 'impulse';
dsp.FIRHalfbandInterpolator и dsp.FIRHalfbandDecimatorThe firhalfband Функция предоставляет несколько другие опции проекта. Однако в большинстве случаев лучше работать непосредственно с dsp.FIRHalfbandInterpolator и dsp.FIRHalfbandDecimator. Эти два системных объектов не только проектируют коэффициенты, но и обеспечивают эффективную реализацию полифазного интерполятора/дециматора. Они поддерживают фильтрацию данных с плавающей точностью двойной/одинарной, а также данных с фиксированной точкой. Они также поддерживают генерацию кода C и HDL, а также оптимизированную генерацию кода ARM ® Cortex ® M и Cortex A.
halfbandInterpolator = dsp.FIRHalfbandInterpolator('SampleRate',Fs, ... 'Specification','Filter order and transition width', ... 'FilterOrder',N,'TransitionWidth',4000); fvtool(halfbandInterpolator,'Fs',2*Fs,'Color','white');
В порядок выполнения интерполяции dsp.FIRHalfbandInterpolator Используется системный объект. Поскольку это многоскоростной фильтр, важно определить, что подразумевается под частотой дискретизации. Для этого и всех других системных объектов частота дискретизации относится к скорости дискретизации входного сигнала. Однако FVTool определяет частоту дискретизации как скорость, с которой работает фильтр. В случае интерполяции вы усиливаете и затем фильтруете (концептуально), поэтому частота дискретизации FVTool должна быть задана как 2 * Fs из-за увеличения дискретизации на 2.
FrameSize = 256; scope = dsp.SpectrumAnalyzer('SampleRate',2*Fs,'SpectralAverages',5); sine1 = dsp.SineWave('Frequency',10e3','SampleRate',Fs, ... 'SamplesPerFrame',FrameSize); sine2 = dsp.SineWave('Frequency',20e3','SampleRate',Fs, ... 'SamplesPerFrame',FrameSize); tic while toc < 10 x = sine1() + sine2() + 0.01.*randn(FrameSize,1); % 96 kHz y = halfbandInterpolator(x); % 192 kHz scope(y); end release(scope);
Заметьте, что спектральные реплики ослаблены приблизительно на 40 дБ, что является примерно ослаблением, обеспечиваемым полуполосе фильтром. Компенсируя групповую задержку в фильтре, можно построить график наложенных входов и интерполированных выборок. Заметьте, что входные выборки остаются неизменными на выходе фильтра. Это происходит потому, что одна из полифазных ветвей полуполоса является чистой ветвью задержки, которая не изменяет входные выборки.
grpDel = 50; n = 0:2:511; stem(n(1:end-grpDel/2),x(1:end-grpDel/2),'k','filled') hold on nu = 0:511; stem(nu(1:end-grpDel),y(grpDel+1:end)) legend('Input samples','Interpolated samples')
В случае десятикратного уменьшения частота дискретизации, заданная в dsp.FIRHalfbandDecimator соответствует скорости дискретизации фильтра, так как вы фильтруете, а затем понижаете частоту (концептуально). Поэтому для дециматоров Fs, заданные в FVTool, не должны умножаться на какой-либо коэффициент.
FrameSize = 256; FsIn = 2*Fs; halfbandDecimator = dsp.FIRHalfbandDecimator('SampleRate',FsIn, ... 'Specification','Filter order and transition width', ... 'FilterOrder',N,'TransitionWidth',4000); fvtool(halfbandDecimator,'Fs',FsIn,'Color','white'); scope = dsp.SpectrumAnalyzer('SampleRate',Fs,'SpectralAverages',5); sine1 = dsp.SineWave('Frequency',10e3','SampleRate',Fs, ... 'SamplesPerFrame',FrameSize); sine2 = dsp.SineWave('Frequency',20e3','SampleRate',Fs, ... 'SamplesPerFrame',FrameSize); tic while toc < 10 x = sine1() + sine2() + 0.01.*randn(FrameSize,1); % 96 kHz y = halfbandInterpolator(x); % 192 kHz xd = halfbandDecimator(y); % 96 kHz scope(xd); end release(scope);
Коэффициенты фильтра могут быть извлечены из интерполятора/дециматора при помощи tf функция.
num = tf(halfbandInterpolator); % Or num = tf(halfbandDecimator);
Вместо определения порядка фильтра и ширины перехода можно спроектировать фильтр минимального порядка, который обеспечивает заданную ширину перехода, а также заданное затухание в полосе задерживания.
Ast = 80; % 80 dB halfbandInterpolator = dsp.FIRHalfbandInterpolator('SampleRate',Fs, ... 'Specification','Transition width and stopband attenuation', ... 'StopbandAttenuation',Ast,'TransitionWidth',4000); fvtool(halfbandInterpolator,'Fs',2*Fs,'Color','white');
Заметьте, что как и во всех интерполяторах, коэффициент усиления полосы пропускания в абсолютных модулях равен коэффициенту интерполяции (2 в случае полуполос). Это соответствует усилению полосы пропускания 6,02 дБ.
Также возможно задать порядок фильтра и затухание в полосе задерживания.
halfbandDecimator = dsp.FIRHalfbandDecimator('SampleRate',Fs, ... 'Specification','Filter order and stopband attenuation', ... 'StopbandAttenuation',Ast,'FilterOrder',N); fvtool(halfbandDecimator,'Fs',Fs,'Color','white');
В отличие от интерполяторов, дециматоры имеют коэффициент усиления 1 (0 дБ) в полосе пропускания.
Полуполосные интерполяторы и дециматоры могут использоваться для эффективной реализации банков фильтров синтеза/анализа. Показанные до настоящего времени полудиапазонные фильтры являются lowpass. С помощью одного дополнительного сумматора можно получить высокочастотный ответ в дополнение к lowpass и использовать два ответа для реализации банка фильтров.
Следующий код моделирует банк квадратурного зеркального фильтра (QMF). Сигнал 8 кГц, состоящий из синусоид 1 кГц и 3 кГц, разделяется на два сигнала 4 кГц с помощью дециматора lowpass/highpass полуполосы. Сигнал lowpass сохраняет синусоиду 1 кГц, в то время как сигнал highpass сохраняет синусоиду 3 кГц (которая сглаживается 1 кГц после понижающей дискретизации). Сигналы затем объединяются назад с набором синтезирующих фильтров с помощью интерполятора с полуполосами. Высокочастотная ветвь преобразует частоту сглаживания 1 кГц синусоиды обратно к 3 кГц. Интерполированный сигнал имеет частоту дискретизации 8 кГц.
Fs1 = 8000; % Units = Hz Spec = 'Filter order and transition width'; Order = 52; TW = 4.1e2; % Units = Hz % Construct FIR Halfband Interpolator halfbandInterpolator = dsp.FIRHalfbandInterpolator( ... 'Specification',Spec, ... 'FilterOrder',Order, ... 'TransitionWidth',TW, ... 'SampleRate',Fs1/2, ... 'FilterBankInputPort',true); % Construct FIR Halfband Decimator halfbandDecimator = dsp.FIRHalfbandDecimator( ... 'Specification',Spec, ... 'FilterOrder',Order, ... 'TransitionWidth',TW, ... 'SampleRate',Fs1); % Input f1 = 1000; f2 = 3000; InputWave = dsp.SineWave('Frequency',[f1,f2],'SampleRate',Fs1, ... 'SamplesPerFrame',1024,'Amplitude',[1 0.25]); % Construct Spectrum Analyzer object to view the input and output scope = dsp.SpectrumAnalyzer('SampleRate',Fs1, ... 'PlotAsTwoSidedSpectrum',false,'ShowLegend',true,'YLimits', ... [-120 30], ... 'Title', ... 'Input Signal and Output Signal of Quadrature Mirror Filter'); scope.ChannelNames = {'Input','Output'}; tic while toc < 10 Input = sum(InputWave(),2); NoisyInput = Input+(10^-5)*randn(1024,1); [Lowpass,Highpass] = halfbandDecimator(NoisyInput); Output = halfbandInterpolator(Lowpass,Highpass); scope([NoisyInput,Output]); end release(scope);
Все представленные до сих пор проекты были оптимальными проектами equiripple. Использование fdesign.interpolator и fdesign.decimator, доступны другие алгоритмы проекта.
Fs = 44.1e3; N = 90; TW = 1000/Fs; % Transition width filtSpecs = fdesign.interpolator(2,'halfband','N,TW',N,TW); equirippleHBFilter = design(filtSpecs,'equiripple','SystemObject',true); leastSquaresHBFilter = design(filtSpecs,'firls','SystemObject',true); kaiserHBFilter = design(filtSpecs,'kaiserwin','SystemObject',true);
Можно сравнить проекты с FVTool. Различные проекты допускают компромиссы между минимальным затуханием в полосе задерживания и большим общим ослаблением.
fvt = fvtool(equirippleHBFilter,leastSquaresHBFilter,kaiserHBFilter, ... 'Fs',2*Fs,'Color','white'); legend(fvt,'Equiripple design','Least-squares design', ... 'Kaiser-window design')
В качестве альтернативы можно задать порядок и затухание в полосе задерживания. Это позволяет проводить компромиссы между общим затуханием в полосе задерживания и шириной перехода.
Ast = 60; % Minimum stopband attenuation filtSpecs = fdesign.interpolator(2,'halfband','N,Ast',N,Ast); equirippleHBFilter = design(filtSpecs,'equiripple','SystemObject',true); kaiserHBFilter = design(filtSpecs,'kaiserwin','SystemObject',true); fvt = fvtool(equirippleHBFilter,kaiserHBFilter,'Fs',2*Fs,'Color','white'); legend(fvt,'Equiripple design','Kaiser-window design')
Конструкции окон Кайзера могут также использоваться в дополнение к проектам equiripple при проектировании фильтра минимального порядка, необходимого для соответствия проектным спецификациям. Фактический порядок для проекта окна Кайзера больше, чем тот, который нужен для проекта equiripple, но общее затухание в полосе задерживания лучше в ответ.
Fs = 44.1e3; TW = 1000/(Fs/2); % Transition width Ast = 60; % 60 dB minimum attenuation in the stopband filtSpecs = fdesign.decimator(2,'halfband','TW,Ast',TW,Ast); equirippleHBFilter = design(filtSpecs,'equiripple','SystemObject',true); kaiserHBFilter = design(filtSpecs,'kaiserwin','SystemObject',true); fvt = fvtool(equirippleHBFilter,kaiserHBFilter,'Fs',Fs,'Color','white'); legend(fvt,'Equiripple Design','Kaiser-window design')
Вместо разработки оконных фильтров Кайзера также возможно получить увеличивающееся затухание в полосе задерживания с модифицированными проектами 'equiripple'.
equirippleHBFilter1 = design(filtSpecs,'equiripple', ... 'StopbandShape','1/f','StopbandDecay',4,'SystemObject',true); equirippleHBFilter2 = design(filtSpecs,'equiripple', ... 'StopbandShape','linear','StopbandDecay',53.333,'SystemObject',true); fvt = fvtool(equirippleHBFilter1,equirippleHBFilter2, ... 'Fs',Fs,'Color','white'); legend(fvt,'Stopband decaying as (1/f)^4','Stopband decaying linearly')
Фильтр highpass halfband может быть получен из фильтра lowpass halfband путем изменения знака каждого второго коэффициента. Кроме того, можно непосредственно спроектировать highpass halfband путем установки свойства 'Type' на 'Highpass'.
filtSpecs = fdesign.decimator(2,'halfband', ... 'Type','Highpass','TW,Ast',TW,Ast); halfbandHPFilter = design(filtSpecs,'equiripple', ... 'StopbandShape','linear','StopbandDecay',53.333,'SystemObject',true); fvt = fvtool(halfbandHPFilter,equirippleHBFilter2,'Fs',Fs,'Color','white'); legend(fvt,'Highpass halfband filter','Lowpass halfband filter')
