功率谱、韦尔奇方法和窗函数答案

【问题标题】：Power spectrum, Welch method and window functions功率谱、韦尔奇方法和窗函数
【发布时间】：2019-11-03 05:32:48
【问题描述】：

我需要使用 Welch 方法找到信号的功率谱。然后使用 Chebyshev、Kaiser 和 Gauss 窗口，找到频率分离方面的最佳窗口。

我写了一个函数来产生带噪声的正弦信号，并对这个信号使用了 Welch 方法，还生成了所有这三个窗口函数。

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

def sig_noise(f,snr,n):
    fs=1000
    dt=1/fs
    tk=n*dt
    A=1
    B=10**(-snr/20)
    t=np.arange(0,tk,dt)
    y=(np.random.rand(n)-0.5)
    z=A*np.sin(2*np.pi*f*t)+B*y
    return z,t

[s,t]=sig_noise(140,-10,64)

f, Pxx_den = sg.welch(s,nperseg=64,noverlap=None) # Welch method, parameters are taken from book

plt.plot(f,Pxx_den)
plt.xlabel("Frequency")
plt.grid()

plt.figure(figsize=(10,10))

chebyshev=sg.chebwin=(64,40) # Chebyshev window

kaiser=sg.kaiser(64,0) # Kaiser window

gauss=sg.gaussian(64,np.std((f,Pxx_den),ddof=0)) # Gauss window

# plots for check

plt.subplot(221)
plt.grid()
plt.title("Gauss")
plt.plot(gauss)
plt.subplot(222)
plt.grid()
plt.title("Chebyshev")
plt.plot(chebyshev)
plt.subplot(223)
plt.grid()
plt.title("Kaiser")
plt.plot(kaiser)

f, Pxx_den = sg.welch(s,nperseg=64,noverlap=None,window=chebyshev)

最后一行代码不起作用，因为我想使用例如 chebyshev 窗口。未知的窗口类型是我得到的答案。谁能告诉我错在哪里？

【问题讨论】：

标签： signal-processing python-3.7 window-functions noise

【解决方案1】：

问题在于以下行中有一个额外的 = 字符：

chebyshev=sg.chebwin=(64,40) # Chebyshev window
#                  ^^^
#                  extra "=" character

这会导致 chebyshev 等于元组 (64,40) 而不是预期的数组（遵循与其他窗口类似的方法）。稍后，当给定一个元组作为 windows 参数时，scipy.signal.welch 尝试通过为元组的第一个元素指定的窗口类型调用 scipy.signal.get_window 来构造窗口。显然64（在您的情况下是元组的第一个元素）不是有效的窗口类型。

因此，解决方法是简单地纠正错字：

chebyshev=sg.chebwin(64,40) # Chebyshev window

或者使用适当构造的元组调用welch：

f, Pxx_den = sg.welch(s,nperseg=64,noverlap=None,window=('chebwin',40)

【讨论】：