matlab fir2 频率响应与幅度响应不对应
matlab fir2 frequency response doesn't correspond to magnitude response
我正在尝试使用 Matlab 中的 fir2
函数设计一个 FIR,并且我正在努力获得与我试图实现的幅度函数相对应的滤波器频率响应。
这是我用于滤波器设计的代码,x
是如图所示的幅度响应
fs = 48000; %sample rate
fny = fs/2; %Nyquist frequency
tabs = 512; %tabs
%frequency resolution 1/3 octave, length=28
f = [0 63 80 100 125 160 200 250 315 400 500 630 800 1000 1250 1600
2000 2500 3150 4000 5000 6300 8000 10000 12500 16000 20000 fny];
fn = f/(fny);
coeffs = fir2(tabs,fn,x);
figure(1);
[freq_response,fc] = freqz(coeffs,tabs);
H = abs(freq_response);
subplot(2,1,1);
semilogx(fc*fny,H, 'r');
hold on;
xlim([20 20000]);
title('frequency response of fir2 filter');
xlabel('frequency [Hz]');
ylabel('magnitude [dB]');
grid on;
set(gca,'XTick',[20 50 100 200 500 1000 2000 5000 10000 20000]);
subplot(2,1,2);
semilogx(f,20*log10(x),'b');
title('desired magnitude function');
xlabel('frequency [Hz]');
ylabel('fir coefficients');
grid on;
xlim([20 20000]);
set(gca,'XTick',[20 50 100 200 500 1000 2000 5000 10000 20000]);
我做错了什么?
首先我们应该观察到函数 freqz
需要滤波器传递函数的分子和分母多项式系数。然后,您还可以指定其他参数,例如频率响应中所需的点数,但这将在 强制性分子和分母多项式系数之后出现。对于没有显式分母的 FIR 滤波器,隐式分母是常数 1
。因此,在您的情况下调用 freqz
的正确方法是:
[freq_response,fc] = freqz(coeffs,1,tabs);
请注意,默认情况下返回的频率是 angular 个频率,这将覆盖 [0,pi]
个范围。要获得以 Hz 为单位的频率,您需要使用 f*fny/pi
重新调整 angular 频率(正如您尝试做的那样,但缺少 pi
因子)。或者,您可以简单地将采样率传递给 freqz
函数:
[freq_response,fc] = freqz(coeffs,1,tabs,fs);
然后您还有一个与图的 y 轴比例相关的问题。鉴于您正在尝试以 dB 为单位绘制幅度响应,您应该将计算出的响应 H
(线性刻度)转换为 dB:
[freq_response,fc] = freqz(coeffs,1,tabs,fs);
H = abs(freq_response);
...
semilogx(fc,20*log10(H), 'r');
在同一张图上绘制期望的和计算的响应应该会得到类似于以下内容的结果:
我正在尝试使用 Matlab 中的 fir2
函数设计一个 FIR,并且我正在努力获得与我试图实现的幅度函数相对应的滤波器频率响应。
这是我用于滤波器设计的代码,x
是如图所示的幅度响应
fs = 48000; %sample rate
fny = fs/2; %Nyquist frequency
tabs = 512; %tabs
%frequency resolution 1/3 octave, length=28
f = [0 63 80 100 125 160 200 250 315 400 500 630 800 1000 1250 1600
2000 2500 3150 4000 5000 6300 8000 10000 12500 16000 20000 fny];
fn = f/(fny);
coeffs = fir2(tabs,fn,x);
figure(1);
[freq_response,fc] = freqz(coeffs,tabs);
H = abs(freq_response);
subplot(2,1,1);
semilogx(fc*fny,H, 'r');
hold on;
xlim([20 20000]);
title('frequency response of fir2 filter');
xlabel('frequency [Hz]');
ylabel('magnitude [dB]');
grid on;
set(gca,'XTick',[20 50 100 200 500 1000 2000 5000 10000 20000]);
subplot(2,1,2);
semilogx(f,20*log10(x),'b');
title('desired magnitude function');
xlabel('frequency [Hz]');
ylabel('fir coefficients');
grid on;
xlim([20 20000]);
set(gca,'XTick',[20 50 100 200 500 1000 2000 5000 10000 20000]);
我做错了什么?
首先我们应该观察到函数 freqz
需要滤波器传递函数的分子和分母多项式系数。然后,您还可以指定其他参数,例如频率响应中所需的点数,但这将在 强制性分子和分母多项式系数之后出现。对于没有显式分母的 FIR 滤波器,隐式分母是常数 1
。因此,在您的情况下调用 freqz
的正确方法是:
[freq_response,fc] = freqz(coeffs,1,tabs);
请注意,默认情况下返回的频率是 angular 个频率,这将覆盖 [0,pi]
个范围。要获得以 Hz 为单位的频率,您需要使用 f*fny/pi
重新调整 angular 频率(正如您尝试做的那样,但缺少 pi
因子)。或者,您可以简单地将采样率传递给 freqz
函数:
[freq_response,fc] = freqz(coeffs,1,tabs,fs);
然后您还有一个与图的 y 轴比例相关的问题。鉴于您正在尝试以 dB 为单位绘制幅度响应,您应该将计算出的响应 H
(线性刻度)转换为 dB:
[freq_response,fc] = freqz(coeffs,1,tabs,fs);
H = abs(freq_response);
...
semilogx(fc,20*log10(H), 'r');
在同一张图上绘制期望的和计算的响应应该会得到类似于以下内容的结果: