时变频率的频谱图与标度图

Question

我正在比较特定信号的 FFT 与 CWT。

我不清楚如何从 CWT 的相应比例图中读取相应的频率和幅度。此外，我的印象是 CWT 相当不精确？

频谱图在预测精确频率方面似乎非常好，但对于 CWT，我尝试了许多不同的小波，结果是一样的。

我监督了什么吗？这不是解决这个问题的合适工具吗？

在下面，您可以找到我的示例源代码和相应的图。

import matplotlib.pyplot as plt
import numpy as np
from numpy import pi as ?
from scipy.signal import spectrogram
import pywt

f_s = 200              # Sampling rate = number of measurements per second in [Hz]
t   = np.arange(-10,10, 1 / f_s) # Time between [-10s,10s].
T1  = np.tanh(t)/2  + 1.0 # Period in [s]
T2  = 0.125               # Period in [s]
f1  = 1 / T1              # Frequency in [Hz]
f2  = 1 / T2              # Frequency in [Hz] 

N = len(t)
x = 13 * np.sin(2 * ? * f1 * t) + 42 * np.sin(2 * ? * f2 * t)

fig, (ax1, ax2, ax3, ax4) = plt.subplots(4,1, sharex = True, figsize = (10,8))

# Signal
ax1.plot(t, x)
ax1.grid(True)
ax1.set_ylabel("$x(t)$")
ax1.set_title("Signal x(t)")

# Frequency change
ax2.plot(t, f1)
ax2.grid(True)
ax2.set_ylabel("$f_1$ in [Hz]")
ax2.set_title("Change of frequency $f_1(t)$")

# Moving fourier transform, i.e. spectrogram
?t = 4 # window length in [s]
Nw = np.int(2**np.round(np.log2(?t * f_s))) # Number of datapoints within window
f, t_, Sxx = spectrogram(x, f_s, window='hanning', nperseg=Nw, noverlap = Nw - 100, detrend=False, scaling='spectrum')
?f  =  f[1] - f[0]
?t_ = t_[1] - t_[0]
t2  = t_ + t[0] - ?t_

im = ax3.pcolormesh(t2, f - ?f/2, np.sqrt(2*Sxx), cmap = "inferno_r")#, alpha = 0.5)
ax3.grid(True)
ax3.set_ylabel("Frequency in [Hz]")
ax3.set_ylim(0, 10)
ax3.set_xlim(np.min(t2),np.max(t2))
ax3.set_title("Spectrogram using FFT and Hanning window")

# Wavelet transform, i.e. scaleogram
cwtmatr, freqs = pywt.cwt(x, np.arange(1, 512), "gaus1", sampling_period = 1 / f_s)
im2 = ax4.pcolormesh(t, freqs, cwtmatr, vmin=0, cmap = "inferno" )  
ax4.set_ylim(0,10)
ax4.set_ylabel("Frequency in [Hz]")
ax4.set_xlabel("Time in [s]")
ax4.set_title("Scaleogram using wavelet GAUS1")

# plt.savefig("./fourplot.pdf")

plt.show()

Answer 1

您的示例波形在任何给定时间仅由几个纯音组成，因此频谱图看起来非常干净且可读。

\n

相比之下，小波图看起来会“混乱”，因为您必须对可能不同尺度（以及频率）的高斯小波求和，以近似原始信号中的每个分量纯音。

\n

短时FFT和小波变换都属于时频变换类别，但具有不同的内核。FFT的内核是纯音，但是您提供的小波变换的内核是高斯小波。FFT 纯音内核将清楚地对应于您在问题中显示的波形类型，但小波内核不会。

\n

我怀疑您对不同小波得到的结果在数值上是否完全相同。它可能只是用肉眼看起来与您绘制的方式相同。看来您将小波用于错误的目的。小波在分析信号方面比绘制信号更有用。小波分析是独特的，因为小波组合中的每个数据点同时编码频率、相位和加窗信息。这允许设计位于时间序列和频率分析之间的连续体的算法，并且可以非常强大。

\n

至于您关于不同小波给出相同结果的说法，这显然是不正确的用于所有小波：我改编了您的代码，它生成了后面的图像。当然，GAUS2 和 MEXH 似乎产生了类似的图（放大到非常近，您会发现它们略有不同），但那是因为二阶高斯小波看起来与墨西哥帽小波相似。

\n

\nimport matplotlib.pyplot as plt\nimport numpy as np\nfrom numpy import pi as \xcf\x80\nfrom scipy.signal import spectrogram, wavelets\nimport pywt\nimport random\n\n\nf_s = 200              # Sampling rate = number of measurements per second in [Hz]\nt   = np.arange(-10,10, 1 / f_s) # Time between [-10s,10s].\nT1  = np.tanh(t)/2  + 1.0 # Period in [s]\nT2  = 0.125               # Period in [s]\nf1  = 1 / T1              # Frequency in [Hz]\nf2  = 1 / T2              # Frequency in [Hz] \n\nN = len(t)\nx = 13 * np.sin(2 * \xcf\x80 * f1 * t) + 42 * np.sin(2 * \xcf\x80 * f2 * t)\n\nfig, (ax1, ax2, ax3, ax4) = plt.subplots(4,1, sharex = True, figsize = (10,8))\n\n\nwvoptions=iter([\'gaus1\',\'gaus2\',\'mexh\',\'morl\'])\n\naxes=[ax1,ax2,ax3,ax4]\n\n\nfor ax in axes:\n    # Wavelet transform, i.e. scaleogram\n    try:\n        choice=next(wvoptions)\n        cwtmatr, freqs = pywt.cwt(x, np.arange(1, 512), choice, sampling_period = 1 / f_s)\n        im = ax.pcolormesh(t, freqs, cwtmatr, vmin=0, cmap = "inferno" )  \n        ax.set_ylim(0,10)\n        ax.set_ylabel("Frequency in [Hz]")\n        ax.set_xlabel("Time in [s]")\n        ax.set_title(f"Scaleogram using wavelet {choice.upper()}")\n    except:\n        pass\n# plt.savefig("./fourplot.pdf")\n\nplt.tight_layout()\nplt.show()\n

\n

\n