在 python 中绘制信号的连续小波变换 (CWT) 的尺度图

信息处理 Python 小波 转换 时频 可视化
2022-02-19 17:53:54

所以我已经学习了几个星期的小波,因为我想在我正在从事的研究项目中使用它们并且我一直在努力掌握它们背后的一般想法。我一直在努力绘制信号 CWT 的尺度图。

如果有人可以通过我的比例图绘图功能来查看我是否正确绘制了这个比例图,我将非常感激。我主要是在努力如何在尺度图中可视化信号的功率水平,或者我现在是否正确地做到了


原始信号和尺度图:

该信号的采样频率为 2048hz,信号长度为 2048 个样本,因此这是我的信号的 1 秒样本。您可以忽略第一个图中的黑线信号。我使用 pywavlets cwt 函数将 cwt 函数应用于原始蓝色信号。

示例信号和尺度图


关于我的 plot_wavelet 函数的问题:

1.在下面的函数中有硬编码的级别值。log2 值位于图右侧的颜色栏中。这些在 ax.contourf 函数中用于创建轮廓线。我应该如何确定要使用多少级别或级别的值?这是我在这里一直在努力解决的主要问题。基本上如何将力量与情节中的颜色联系起来。

2.在我的尺度图的这一点上,我假设接近值 0 的区域具有最大的功率。当我将位置与我在原始信号图中看到的位置进行比较时,这是有道理的。这是一个正确的假设吗?我想知道如何使图形的分辨率更好,或者这甚至可能吗?如果我使用更长范围的刻度,我想我可以提高频带分辨率?这是不是正确的思路。

3.我应该如何在比例图中显示我的 y 轴?pywt.cwt 函数返回频率,但我猜它们更像是频率范围?我在理解如何准确解释 y 轴方面有点挣扎......我已经阅读了很多不同的论文/教程,我想我对如何关联比例和频率或如何以这种方式绘制它。

plot_wavelet 函数代码:

def plot_wavelet(ax, time2, signal, scales, waveletname = 'cmor', 
                 cmap =plt.cm.seismic, title = '', ylabel = '', xlabel = ''):
    dt=time2
    coefficients, frequencies = pywt.cwt(signal, scales, waveletname, dt)

    power = (abs(coefficients)) ** 2
    period = frequencies
    levels = [0.015625,0.03125,0.0625, 0.125, 0.25, 0.5, 1]
    contourlevels = np.log2(levels) #original
    time=range(2048)

    im = ax.contourf(time, np.log2(period), np.log2(power), contourlevels, extend='both',cmap=cmap)


    ax.set_title(title, fontsize=20)
    ax.set_ylabel(ylabel, fontsize=18)
    ax.set_xlabel(xlabel, fontsize=18)
    yticks = 2**np.arange(np.ceil(np.log2(period.min())), np.ceil(np.log2(period.max())))    
    ax.set_yticks(np.log2(yticks)) #original
    ax.set_yticklabels(yticks) #original
    ax.invert_yaxis()
    ylim = ax.get_ylim()

    cbar_ax = fig.add_axes([0.95, 0.5, 0.03, 0.25])
    fig.colorbar(im, cax=cbar_ax, orientation="vertical")

    return yticks, ylim


创建上述两个图的代码,请记住您可以忽略第一个图中的黑色信号。我对原始蓝色信号执行了 CWT。


xrange=list(range(2048))
fig, ax = plt.subplots(figsize=(12,8))
ax.plot(xrange,signal, color="b", alpha=0.5, label='original signal')
#rec = lowpassfilter(signal, 0.4)
ax.plot(xrange,rec, 'k', label='DWT smoothing}', linewidth=2)
ax.legend()
ax.set_title('Removing High Frequency Noise with DWT', fontsize=18)
ax.set_ylabel('Signal Amplitude', fontsize=16)
ax.set_xlabel('Sample No', fontsize=16)
plt.margins(0)
plt.show()


scale_range = np.arange(2, 50) # number of scales
fig, ax = plt.subplots(figsize=(12, 8))
plot_wavelet(ax=ax, time2=sp, signal=signal, scales=scale_range,waveletname='cmor1.5-1.0',
             title = "CWT of Signal", ylabel = ylabel, xlabel = xlabel)
plt.show()


参考链接:

链接到我从哪里获得我正在关注的示例以及绘图代码的来源

链接到 pywavelets cwt 函数

0个回答
没有发现任何回复~
其它你可能感兴趣的问题