小波变换原理简要实现和分析

傅里叶变换是频谱分析的重要工具，适用于周期性平稳信号，但是对于非平稳信号而言，效果较差。

其实质是卷积求正弦波相关性，由于正弦波不是能量有限信号，对无差别做一个全部卷积，算出来的结果不考虑时间。

小波变换也是频谱分析的重要工具，基函数为能量有限信号，也可以进行时域分析。针对不同的信号，可以更换不同的基函数。

建议看看：形象易懂讲解算法I——小波变换 - Alexander - 博客园 (cnblogs.com)

一般有母小波、父小波。母小波就是基函数，父小波是尺度函数，他们都是正交归一化的基向量。父小波是信号在不同尺度空间下的近似。如下所示，为母小波的系数求法，a是尺度，尺度与频率是反比关系，为了保证归一化，要除以根号a。

关于父小波：小波变换完美通俗讲解系列之（二） - 知乎 (zhihu.com)

上式的第一项是父小波，第二项则是母小波。

下面是cwt的一个小小实现，抄自(26条消息) python连续小波分析CWT_python 小波变换_赵孝正的博客-CSDN博客

(26条消息) 脑电信号处理学习笔记（一）——pywt_pywt.wavedec_牛杂师傅的博客-CSDN博客

import numpy as np
import pywt
import matplotlib
matplotlib.use('TkAgg')
import matplotlib.pylab as plt

sampling_rate = 1024  # 采样频率
t = np.arange(0, 1.0, 1.0/sampling_rate)  # 0-1.0之间的数，步长为1.0/sampling_rate
f1 = 100  # 频率
f2 = 200
f3 = 300
data = np.piecewise(t, [t < 1, t < 0.8, t < 0.3],
                    [lambda t: np.sin(2 * np.pi * f1 * t),
                     lambda t: np.sin(2 * np.pi * f2 * t),
                     lambda t: np.sin(2 * np.pi * f3 * t)]
                    )
wavename = "cgau8"  # 小波函数
totalscal = 256     # totalscal是对信号进行小波变换时所用尺度序列的长度(通常需要预先设定好)
fc = pywt.central_frequency(wavename)  # 计算小波函数的中心频率
cparam = 2 * fc * totalscal  # 常数c
scales = cparam/np.arange(totalscal, 1, -1)  # 为使转换后的频率序列是一等差序列，尺度序列必须取为这一形式（也即小波尺度）
[cwtmatr, frequencies] = pywt.cwt(data, scales, wavename, 1.0/sampling_rate)  # 连续小波变换模块

plt.figure(figsize=(8, 4))
plt.subplot(211)  # 第一整行
plt.plot(t, data)
plt.xlabel(u"time(s)")
plt.title(u"300Hz 200Hz 100Hz Time spectrum")
plt.subplot(212)  # 第二整行

plt.contourf(t, frequencies, abs(cwtmatr))
plt.ylabel(u"freq(Hz)")
plt.xlabel(u"time(s)")
plt.subplots_adjust(hspace=0.4) #调整边距和子图的间距 hspace为子图之间的空间保留的高度，平均轴高度的一部分
plt.show()

关于如何选取基函数，请移步到：(26条消息) 小波滤波小波基、阈值、阈值函数的选择_小波基函数的选取_狂奔的蜗牛的博客-CSDN博客

小波基的选取 - 知乎 (zhihu.com)

此外提供一些方法查看小波形态

import matplotlib.pyplot as plt
wavelet=pywt.ContinuousWavelet('gaus8')
[psi,xval] = wavelet.wavefun( level=10)#或者写length=1024
plt.plot(xval,psi)
plt.title("Gaussian Wavelet of order 1024")
plt.show()