在 Python 中使用双线性插值连续到离散答案

【问题标题】：Continuous to discrete using bilinear interpolation in Python在 Python 中使用双线性插值连续到离散
【发布时间】：2018-12-05 16:01:54
【问题描述】：

我在频域中有一个一阶低通滤波器 (LPF)，我想将它数字化。我正在比较频率响应图以进行测试，但我得到了奇怪的结果...

虽然非常基本，但我无法通过阅读 scipy.signal.bilinear 帮助页面或网络来正确理解。

num,den：S平面的分子和分母
b, a: 我希望它是 b, a 数字差分方程滤波器 (IIR) 的系数，形状为：Y[n] = a0 *X[n] + a1*X[N-1] + ... - b1*Y[n-1] ...

这是一个代码示例：

Fs = 48000.0
f = 2 * np.logspace(1,4,1024)

num =  [0 , 1]
den = [0.001 , 1]

tmp, H = sig.freqs(num, den, worN=1024)
b, a = sig.bilinear(num, den, 1.0)
tmp, Hd = sig.freqz(b,a, worN=1024)

plt.semilogx(f, 20*np.log10(np.abs(H)))
plt.semilogx(f, 20*np.log10(np.abs(Hd)))

我做错了什么？

【问题讨论】：

标签： python signal-processing

【解决方案1】：

问题是您在绘图时没有在 x 轴上使用 tmp，而且 freqz 以弧度/样本给出了归一化的 tmp 向量：

import numpy as np
import scipy.signal as sig
import matplotlib.pyplot as plt

Fs = 48000
num =  [0 , 1000]
den = [1 , 1000]

w1, H = sig.freqs(num, den, worN=1024)
b, a = sig.bilinear(num, den, Fs)
w2, Hd = sig.freqz(b, a, worN=1024)

fig = plt.figure()
plt.title('Filter frequency response')
plt.semilogx(w1, 20*np.log10(np.abs(H)),'b')
plt.semilogx(w2*Fs, 20*np.log10(np.abs(Hd)),'k')
plt.ylabel('magnitude [dB]')
plt.xlabel('frequency [Hz]')
plt.grid()
plt.axis('tight')
plt.xlim([0.001, Fs/2])
plt.show()

这段代码运行良好。希望对您有所帮助。

【讨论】：

【解决方案2】：

是的，这很有帮助。然后我将其带到下一步，以便两个“w”向量的长度相同并覆盖所有有趣的采样点，如下所示：

Fs = 48000.0
f = 2 * np.logspace(1,4,1024)
w = 2 * np.pi * f

Snum =  [0 , 1]
Sden = [0.001 , 1]

w1, H = sig.freqs(Snum, Sden, worN=w)
b, a = sig.bilinear(Snum, Sden, Fs)
w2, Hd = sig.freqz(b,a, worN=w1/Fs)

另外，如果呈现的图表以 Hz 为单位（不是 Rad/Sec0），w 应除以 2*PI：

plt.semilogx(w1/np.pi/2, 20*np.log10(np.abs(H)), 'r')

【讨论】：