时间戳

首页 > 解决方案 > 在 Python 中绘制 DWT 比例图

python - 在 Python 中绘制 DWT 比例图

问题描述

我有一个来自磁性探测器的信号，我有兴趣分析它，我使用 wavedec() 进行了信号分解

coeffs = pywt.wavedec(dane_K180_40['CH1[uV]'], 'coif5', level=5)

我收到的分解系数如下：

cA1, cD5, cD4, cD3, cD2, cD1 = coeffs

这些是具有各种长度的 ndarrays 对象。cD1 是 (1519,) cD2 是 (774,) 等等。不同长度的数组是我的主要障碍。 系数

我的问题：

我必须制作DWT Scaleogram，我不能强调我已经尽了最大努力却做不到。

最好的方法是什么？使用 matpllotlib 的imshow()如下：

plt.imshow(np.abs([cD5, cD4, cD3, cD2, cD1]), cmap='bone', interpolation='none', aspect='auto')

给我一个错误

TypeError: dtype 对象的图像数据不能转换为浮点数

因为我不是 python 专家，所以我试图用谷歌搜索它，并且我试图将 ndarrays 更改为浮动。

绘制 scaleogram、matshow、pcolormesh 的最佳方法是什么？;D

标签： pythonsignal-processingwavelet-transformpywavelets

解决方案

基本上，每个 cDi 数组的样本量是前一个数组的一半（每个母小波都不是这种情况！），所以我创建了一个 2D numpy 数组，其中第一个元素是“全部”样本量，对于每个后续级别我都会重复采样2^level次数，以便最终结果是一个矩形块。您可以选择是否要将 Y 轴绘制为线性或对数刻度。

# Create signal
xc = np.linspace(0, t_n, num=N)
xd = np.linspace(0, t_n, num=32)
sig = np.sin(2*np.pi * 64 * xc[:32]) * (1 - xd)
composite_signal3 = np.concatenate([np.zeros(32), sig[:32], np.zeros(N-32-32)])

# Use the Daubechies wavelet
w = pywt.Wavelet('db1')
# Perform Wavelet transform up to log2(N) levels
lvls = ceil(log2(N))
coeffs = pywt.wavedec(composite_signal3, w, level=lvls)
# Each level of the WT will split the frequency band in two and apply a
# WT on the highest band. The lower band then gets split into two again,
# and a WT is applied on the higher band of that split. This repeats
# 'lvls' times.
#
# Since the amount of samples in each step decreases, we need to make
# sure that we repeat the samples 2^i times where i is the level so
# that at each level, we have the same amount of transformed samples
# as in the first level. This is only necessary because of plotting.
cc = np.abs(np.array([coeffs[-1]]))
for i in range(lvls - 1):
    cc = np.concatenate(np.abs([cc, np.array([np.repeat(coeffs[lvls - 1 - i], pow(2, i + 1))])]))

plt.figure()
plt.xlabel('Time (s)')
plt.ylabel('Frequency (Hz)')
plt.title('Discrete Wavelet Transform')
# X-axis has a linear scale (time)
x = np.linspace(start=0, stop=1, num=N//2)
# Y-axis has a logarithmic scale (frequency)
y = np.logspace(start=lvls-1, stop=0, num=lvls, base=2)
X, Y = np.meshgrid(x, y)
plt.pcolormesh(X, Y, cc)

use_log_scale = False

if use_log_scale:
    plt.yscale('log')
else:
    yticks = [pow(2, i) for i in range(lvls)]
    plt.yticks(yticks)

plt.tight_layout()
plt.show()

推荐阅读