matlab与python scipy.signal中的freqs,freqz频率分析函数,输出的w,有时候是角频率,有时候是真实频率,容易搞混,这里对比一下。
1) freqs:
matlab, 角频率,rad.s
python, 角频率 rad/s ,只能是角频率。
2) freqz
matlab, 形式为 [h,w] = freqz(b,a,n) 角频率
形式为 [h,f] = freqz(___,n,fs) 时,频率输出形式f为Hz形式,fs为采样频率
python scipy 中 w,h =freqz(b,a,worN,fs) , w的单位与输入fs相同,fs为归一化角频率时,w也为角频率,fs为采样频率,单位Hz时,w也为Hz。
3) 角频率范围的区别:
freqs中的角频率是现实中的量,可以很大,比如1000Hz,对应的角频率为1000*2*pi ; freqz中的角频率是数字化的,一般使用时是归一化的,范围在 0,2*pi之间。
4) 角频率与Hz频率转化:
freqs的w结果要想用Hz,显示,可以先 w/2/pi 转化为 Hz 频率; freqz中的角频率如果要转化为具体的频率, 因为他是归一化的,用 0~ pi 的范围代表 0- fs/2 的频率范围,可以用 f=( w/pi)*(fs/2) 转化为Hz频率
freqs 是角频率w的单位 rad/s,想要变成Hz, 显示时使用 f = w/2/pi
模拟的freqs不存在归一化。
a = [1 0.4 1];
b = [0.2 0.3 1];
w = logspace(-1,1);
h = freqs(b,a,w);
mag = abs(h);
phase = angle(h);
phasedeg = phase*180/pi;
subplot(2,1,1)
loglog(w,mag)
grid on
xlabel('Frequency (rad/s)')
ylabel('Magnitude')
subplot(2,1,2)
semilogx(w,phasedeg)
grid on
xlabel('Frequency (rad/s)')
ylabel('Phase (degrees)')
freqs 输出的w也是rad/s,也只能是rad/s 角频率。不过这个不是归一化的。模拟的freqs不存在归一化。
w : ndarray
The angular frequencies at which `h` was computed.
b = [1]
a = [0.125 ,1]
from scipy.signal import bilinear,freqs,freqz
import matplotlib.pyplot as plt
import numpy as np
wf=np.logspace(-1, 4, 1000)
w,h = freqs(b,a,wf)
plt.semilogx(w,20*np.log10(np.abs(h)))
plt.xlabel('rad/s')
如何转化为Hz显示,就是x坐标轴 除以 2*pi
b = [1]
a = [0.125 ,1]
from scipy.signal import bilinear,freqs,freqz
import matplotlib.pyplot as plt
import numpy as np
wf=np.logspace(-1, 4, 1000)
w,h = freqs(b,a,wf)
plt.semilogx(w/2/np.pi,20*np.log10(np.abs(h)))
plt.xlabel('Hz')
1) 函数形式为
[h,w] = freqz(b,a,n)
时,w输出为角频率,且归一化,即最大的角频率为 pi (对应fs/2,归一化处理) 。(n为输出的点的个数)和freqs中
2) 函数形式为
[h,f] = freqz(___,n,fs)
时,频率输出形式f为Hz形式,fs为采样频率。
b=1;
a=[0.125 1];
fs=2000;
[bz,az] = bilinear(b,a,fs);
[hz0,wz0] = freqz(bz,az,100); % 100是 n,输出点的个数
fz0= wz0/pi*fs/2; % 将 归一化的rad/s 转化为 实际的采样频率
plot(fz0,20*log10(abs(hz0)),'o');
xlabel('Hz')
hold on
[hz1,fz1] = freqz(bz,az,100,fs);
plot(fz1,20*log10(abs(hz1)),'r-','linewidth',3);
legend('wz0' , 'fz1')
hold off
结果:
freqz(b,a=1, worN=512, whole=False, plot=None, fs=2*pi, include_nyquist = False,)
Returns
-------
w : ndarray
The frequencies at which `h` was computed, in the same units as `fs`.
By default, `w` is normalized to the range [0, pi) (radians/sample).
w的单位和输入fs的单位相同,如果fs是用的 rad/s则返回w也是rad/s, 若输入fs的单位是 Hz,那么输出的w单位也是Hz。
代码部分
from scipy.signal import bilinear,freqs,freqz
import matplotlib.pyplot as plt
import numpy as np
b = [1]
a = [0.125 ,1]
wf=np.logspace(-1,4,1000 )
w,h = freqs(b,a,wf)
plt.semilogx(w/2/np.pi,20*np.log10(np.abs(h)),'o',label='freqs')
plt.xlabel('Hz')
plt.ylabel('dB')
fs=5000
bz,az = bilinear(b,a,fs)
worN=np.logspace(-1,4,2000)
idx_end = np.nonzero(worN<=fs/2)[0][-1]
z = freqz(bz,az,worN=worN[0:idx_end],fs=fs)
plt.semilogx(z[0],20*np.log10(z[1]),'-',label='freqz')
plt.legend()
手机扫一扫
移动阅读更方便
你可能感兴趣的文章