【问题标题】:Frequency Response Scipy.signal频率响应 Scipy.signal
【发布时间】:2016-10-16 11:31:42
【问题描述】:

我正在学习数字信号处理以实现过滤器,并正在使用 python 轻松实现测试想法。所以我才开始使用 scipy.signal 库来查找不同滤波器的脉冲响应和频率响应。

目前,我正在阅读 Paul A. Lynn (1992) 所著的“数字信号、处理器和噪声”一书(并发现它是学习这些东西的绝佳资源)。在本书中,他们有一个过滤器,其传递函数如下所示:

我将分子和分母除以,得到以下等式:

然后我使用 Scipy 实现了这一点:

NumeratorZcoefs = [1, -1, 1, -1]
DenominatorZcoefs = [1, 0.54048, -0.62519, -0.66354, 0.60317, 0.69341]

FreqResponse = scipy.signal.freqz(NumeratorZcoefs, DenominatorZcoefs)
fig = plt.figure(figsize = [8, 6])
ax = fig.add_subplot(111)
ax.plot(FreqResponse[0], abs(np.array(FreqResponse[1])))
ax.set_xlim(0, 2*np.pi)
ax.set_xlabel("$\Omega$")

并生成如下图:

但是在书中显示的频率响应如下:

它们的形状相同,但 2.3 和 0.5 处的峰值比率对于 2 个图来说非常不同,有人可以建议这是为什么吗?

编辑:

除此之外,我刚刚实现了一个手动计算频率响应的函数(通过计算函数极点和零点的距离),我得到了与 scipy.signal 生成的图相似的比率,但是数字不一样,有人知道为什么会这样吗?

实现如下:

def H(omega):
    z1 = np.array([0,0]) # zero at 0, 0
    z2 = np.array([0,0]) # Another zero at 0, 0
    z3 = np.array([0, 1]) # zero at i
    z4 = np.array([0, -1]) # zero at -i
    z5 = np.array([1, 0]) # zero at 1

    z = np.array([z1, z2, z3, z4, z5])

    p1 = np.array([-0.8, 0])
    p = cmath.rect(0.98, np.pi/4)
    p2 = np.array([p.real, p.imag])
    p = cmath.rect(0.98, -np.pi/4)
    p3 = np.array([p.real, p.imag])
    p = cmath.rect(0.95, 5*np.pi/6)
    p4 = np.array([p.real, p.imag])
    p = cmath.rect(0.95, -5*np.pi/6)
    p5 = np.array([p.real, p.imag])

    p = np.array([p1, p2, p3, p4, p5])

    a = cmath.rect(1,omega)
    a_2dvector = np.array([a.real, a.imag])

    dz = z-a_2dvector
    dp = p-a_2dvector

    dzmag = []
    for dis in dz:
           dzmag.append(np.sqrt(dis.dot(dis)))

    dpmag = []
    for dis in dp:
           dpmag.append(np.sqrt(dis.dot(dis)))        

    return(np.product(dzmag)/np.product(dpmag))

然后我像这样绘制频率响应:

omegalist = np.linspace(0,2*np.pi,5000)
Hlist = []

for omega in omegalist:
    Hlist.append(H(omega))

fig = plt.figure()
ax = fig.add_subplot(111)
ax.plot(omegalist, Hlist)
ax.set_xlabel("$\Omega$")
ax.set_ylabel("$|H(\Omega)|$")

并得到以下情节:

【问题讨论】:

    标签: python scipy signal-processing


    【解决方案1】:

    SciPy 生成的频率响应是正确的。无论如何,我不会相信这本书的人物似乎是手绘的。

    如果您想“手动”找到频率响应,这可以通过定义一个返回原始 Z 变换并在单位圆上对其进行评估的函数来简单地完成,如下所示

    def H(z):
        num = z**5 - z**4 + z**3 - z**2
        denom = z**5 + 0.54048*z**4 - 0.62519*z**3 - 0.66354*z**2 + 0.60317*z + 0.69341
        return num/denom
    
    import numpy as np
    import matplotlib.pyplot as plt
    
    w_range = np.linspace(0, 2*np.pi, 1000)
    plt.plot(w_range, np.abs(H(np.exp(1j*w_range))))
    

    结果和 SciPy 完全一样。

    【讨论】:

      猜你喜欢
      • 2013-06-06
      • 1970-01-01
      • 1970-01-01
      • 1970-01-01
      • 1970-01-01
      • 1970-01-01
      • 2011-12-02
      • 1970-01-01
      • 2011-03-11
      相关资源
      最近更新 更多