企业官网建设流程全解析

简介：面向小波变换入门者的MATLAB连续小波变换快捷程序，适合正在进行信号处理或图像处理课程学习的学生，也适合刚接触小波理论的科研人员快速验证想法。程序聚焦连续小波变换的核心实现，压缩包仅由1个.m脚本组成，大小仅795B，结构单纯便于逐行研读，省去了复杂工程配置的干扰。截至目前已有151人学习，说明其在同类入门资源中受到一定关注。下载并运行该脚本后，读者能直观看到连续小波变换的编程步骤，包括小波基选择、尺度与平移参数设置、系数计算与结果展示，有助于将教材公式转化为可执行的代码认知。无论是快速验证信号分解效果，还是对比不同小波基函数的差异，这个脚本都能提供直观的实验载体，可作为课程报告、结课设计或课题探索的基线程序。

1. 从 wavelet.rar 到连续小波变换：先搞懂它在算什么

如果你刚下载了一个叫 wavelet.rar 的压缩包，解压后大概率会看到一堆 .m 或 .py 文件、几张像频谱图的 png，以及一个写着“连续小波变换”的 README。最容易卡住的不是调库，而是不知道程序里那个三维图到底在算什么。

连续小波变换（CWT）把一维信号投影到一组由小波基伸缩平移得到的函数上，输出随尺度 a 和时间 b 变化的系数矩阵。它和短时傅里叶变换最大的区别是窗口自动随频率变化：低频窗宽、高频窗窄，正好匹配非平稳信号对时间分辨率的需求。适合做振动故障诊断、脑电分析、语音处理，也常被用来给图像增强做多尺度分解打基础。

下面沿着“公式→代码→调参→排错”的顺序，把连续小波变换从能查到定义到能跑出靠谱谱图的过程讲清楚。回头再看 wavelet.rar 里的代码就不晕了。

2. 连续小波变换的公式拆解与变尺度窗口的必然性

2.1 CWT 公式里每个符号在代码里对应什么

连续小波变换的常见写法是：

$$X_w(a,b)=\frac{1}{\sqrt{|a|}}\int_{-\infty}^{+\infty} x(t)\psi^*\left(\frac{t-b}{a}\right)dt$$

这里的 $a$ 叫尺度，$b$ 叫平移，$\psi(t)$ 是小波基函数，$*$ 表示复共轭。在 PyWavelets 里，一行pywt.cwt(data, scales, wavelet, sampling_period)几乎把这个公式原样搬进了参数：data对应 $x(t)$，scales数组里的每个值就是一个 $a$，wavelet对应 $\psi(t)$，sampling_period负责把离散采样点索引换算成物理时间。返回的coefs[i][j]就是第 $i$ 个尺度、第 $j$ 个采样时刻的复系数。

看公式时先抓住三个点：尺度小对应高频，小波基被压缩；尺度大对应低频，小波基被拉伸。$b$ 只负责平移，不改变小波基形状。$\psi(t)$ 的均值必须为零，否则积分在无穷远处不收敛，这被称为容许条件（admissibility condition）。Morlet、墨西哥帽、高斯小波都满足这个条件，这也是为什么它们能直接用作 CWT 小波基。

最小验证代码：

import numpy as np import pywt fs = 1000.0 t = np.linspace(0, 1, 1000, endpoint=False) x = np.sin(2 * np.pi * 50 * t) + 0.5 * np.sin(2 * np.pi * 120 * t) scales = np.arange(1, 100, 1) coefs, freqs = pywt.cwt(x, scales, 'morl', sampling_period=1/fs) print(coefs.shape) # (99, 1000) print(freqs[:3]) # 最高频的3个频率值

coefs.shape的第一维等于len(scales)，第二维等于信号长度。freqs的长度也等于len(scales)，但顺序是高频到低频，所以freqs[0]对应第一个尺度，也就是最小时尺度。这里故意用np.arange(1,100,1)做演示，实际工程里应该用对数间隔的尺度序列，后面第 3 章会说明。

2.2 为什么连续小波变换不是“换了个窗”的短时傅里叶

短时傅里叶变换（STFT）用固定长度窗函数截信号。窗长定了，时间分辨率和频率分辨率的乘积就受海森堡不确定性原理限制，想提高时间分辨率就必须牺牲频率分辨率，反之亦然。连续小波变换用 $a$ 去拉伸或压缩 $\psi(t)$：分析高频时等效窗很短，能定位突变的时刻；分析低频时等效窗很长，能区分接近的频率成分。这种变尺度窗口正好匹配语音、振动、生物电信号里常见的高频瞬态叠加低频趋势的结构。

图像增强里反复提到的“多尺度”也来自这个特性。图像的边缘、纹理在不同尺度下呈现不同粗细，单一固定窗滤波很难同时保留细纹理和粗轮廓。小波族天然把尺度维保留下来，后续增强可以按尺度分别处理。用 Python 做图像增强时，工程上更常用二维离散小波变换（wavedec2），因为它计算快、可逆、冗余低，但选小波基、控制分解层数的逻辑与 CWT 的尺度选择完全同源。

2.3 拿到 wavelet.rar 后先别跑主程序，先做三件事

我拿到这类资料包，不会直接双击运行主脚本。里面常有数据，也常有画图代码，但参数往往是写死的。常见做法是先做三件事：确认信号采样率；确认数据是否包含 NaN 或直流漂移；用小波变换自带的自测函数验证频率轴。

这里给一个简单的自测函数：

def inspect_signal(filename, fs): """读入信号并检查基本质量，返回去均值后的数据""" x = np.loadtxt(filename, delimiter=',') # 具体读取按文件格式调整 x = np.nan_to_num(x) x = x - np.mean(x) # CWT 对直流不敏感，去均值可以减少边界伪影 print(f"length={len(x)}, fs={fs}, duration={len(x)/fs:.2f}s") return x

注意：np.loadtxt只适用纯数值文本格式。如果wavelet.rar里是.mat文件，就要用scipy.io.loadmat；如果是.xlsx，用pandas.read_excel。不要盲目照抄读取代码。去均值这步很关键，因为连续小波变换本质是通过积分提取信号能量，直流分量会造成尺度轴两端出现明显的水平亮带，干扰频率峰值的判断。

3. 用 Python 和 PyWavelets 跑通连续小波变换的最小闭环

3.1 最小可运行版本：一条信号、一组尺度、一张图

在确认数据质量后，最小闭环包含四个动作：构造测试信号、选尺度序列、调用pywt.cwt、画出模值谱图。下面是完整代码：

import numpy as np import pywt import matplotlib.pyplot as plt fs = 1000.0 t = np.arange(0, 2, 1/fs) # 频率从 100 Hz 开始衰减的振荡信号，模拟冲击响应 x = np.sin(2 * np.pi * 100 * t) * np.exp(-2 * t) scales = np.geomspace(5, 200, 64) # 对数尺度 coefs, freqs = pywt.cwt(x, scales, 'cmor1.5-1.0', sampling_period=1/fs) plt.figure(figsize=(8, 5)) plt.pcolormesh(t, freqs, np.abs(coefs), shading='auto', cmap='jet') plt.yscale('log') plt.xlabel('Time (s)') plt.ylabel('Frequency (Hz)') plt.colorbar(label='Magnitude') plt.show()

代码里的np.geomspace(5, 200, 64)是尺度范围的工程默认写法。尺度序列如果太密，相邻频率差别不明显，图像上会产生横向过度平滑；太疏又会漏掉窄带信号。64 个尺度用于预览，正式分析用 128。plt.pcolormesh会自动处理不均匀的频率轴，比imshow少踩坐标顺序的坑。

3.2 scales 和 freqs 的换算：写给不懂伪频率的人

pywt.cwt在没传sampling_period时，freqs返回的不是物理频率，而是以采样间隔为单位的“伪频率”。很多教程省略这个参数，导致画出来的频率轴数值与信号实际频率对不上。正确的换算关系是：

物理频率 = 小波中心频率 / (尺度 × 采样周期)

可以用pywt.scale2frequency验证：

fc = 1.0 # cmor1.5-1.0 中的中心频率 freq_expected = pywt.scale2frequency('cmor1.5-1.0', scales) / (1/fs) print(np.max(np.abs(freqs - freq_expected)))

这个差值正常情况下小于1e-12。如果你用的是morl，它的中心频率是 0.8125 而不是 1.0，没注意这点就会把所有频率成比例算偏。所以做定量分析时，优先选择中心频率明确的小波，比如cmorB-F格式，其中B是带宽，F是中心频率。

3.3 “小波变换图像增强python”里为什么更常看到 wavedec2

连续小波变换对一维信号的时频分析效果直观，但二维图像直接做 CWT 计算量非常大，而且连续小波基之间有高冗余。实际搜“小波变换图像增强python”时，绝大多数可运行的代码使用的是二维离散小波变换 DWT。它本质上是把 CWT 的尺度按 2 的幂次离散化，用滤波器组实现快速分解与重建。

import pywt import numpy as np from PIL import Image img = np.array(Image.open('lena.png').convert('L')).astype(float) coeffs = pywt.wavedec2(img, 'db4', level=3) cA3, (cH3, cV3, cD3), (cH2, cV2, cD2), (cH1, cV1, cD1) = coeffs threshold = 0.1 * np.max(np.abs(cH1)) for cH in (cH3, cH2, cH1): cH[cH] = np.sign(cH[cH]) * np.maximum(np.abs(cH[cH]) - threshold, 0) img_enhanced = pywt.waverec2(coeffs, 'db4') img_enhanced = np.clip(img_enhanced, 0, 255).astype(np.uint8)

说明：wavedec2里level=3对应三个尺度。cA3是低频近似，cH/cV/cD分别是水平、垂直、对角高频细节。这里的阈值门限用的是全局阈值 0.1 倍最大值，实际项目中要按不同层单独设置，因为前面尺度的噪声能量分布不一样。db4的消失矩是 4，对图像边缘有较好的稀疏表示。如果你手里wavelet.rar里的 CWT 代码无法直接处理图像，可以先对图像的每一行做 CWT，再按行拼接出二维时频图，但更标准的路径还是切到 DWT。

3.4 什么时候坚持 CWT，什么时候换成 DWT

CWT 适合分析长度几百到几百万点的一维信号；图像这类二维数据，即使只做边缘增强，也建议先 DWT 粗分解。这样代码性能和结果可复现性都比强行 CWT 好。

4. 连续小波变换的 5 个必调参数与实际选参顺序

4.1 小波基的参数格式：从morl到cmorB-F

pywt.cwt的wavelet参数可以直接传字符串，比如'morl'、'cmor1.5-1.0'、'gaus8'。字符串里的数字是有含义的：cmor1.5-1.0表示复数 Morlet 小波，带宽 1.5，中心频率 1.0；gaus8表示高斯小波的 8 阶导数。选型没有绝对标准，但常见经验是：

带宽 B 和中心频率 F 一旦变化，同一个尺度数组对应的物理频率就变了。设置cmor时，带宽太小时域支撑长，边界效应会占掉很大面积；带宽太大则频率分辨率粗糙。我通常从cmor1.5-1.0开始，再根据谱图效果微调。

4.2 尺度范围的计算：目标是覆盖目标频段

尺度范围和采样率、小波中心频率的关系：

scales = fc × fs / f，其中 fs 是采样率，f 是希望分析的物理频率。

所以想分析 5 Hz 到 200 Hz、采样率 1000 Hz、中心频率 1.0 的信号，尺度下限是 1×1000/200=5，尺度上限是 1×1000/5=200。这个简单换算能避免盲目把尺度取到 1000 以上。用代码落地：

fc = 1.0 fs = 1000.0 fmin, fmax = 5.0, 200.0 scale_min = fc * fs / fmax # 5.0 scale_max = fc * fs / fmin # 200.0 scales = np.geomspace(scale_min, scale_max, 128) print(scale_min, scale_max)

注意：fmax不要超过奈奎斯特频率fs/2。超过之后，采样数据里已经没有有效能量，只会增加高频边界伪影。如果信号有 50 Hz 工频干扰，可以按需要把频段下限定在 60 Hz 以上，而不是把宽带噪声一起纳入。

4.3 尺度点数：不是越密越好

尺度序列点数是 CWT 唯一的“分辨率旋钮”。点数越多，coefs矩阵行数越多，计算时间和内存线性上升。对 1 秒采样率 1000 的信号，128 个尺度运行时间通常在几十毫秒到几百毫秒之间；取 512 个尺度图像会显得平滑，但不会产生新信息。

经验值：预览 64，验证参数 128，做峰值搜索 256。如果信号本身只有一个窄带频率，64 点足够；如果存在两个频率接近的振荡，尺度点数要足以让它们在频率轴上分开，此时可以用解析信号来检验分辨率，而不是盲目增加点数。

4.4 边界效应和去趋势：最容易让谱图骗人的两个点

CWT 在信号两端会引入边界伪影，因为小波支撑可能超出数据范围。常见的缓解方法是先对信号做对称延拓（mode='periodization'或'symmetric'），或者分析时把两端各截掉最小尺度对应的小波支撑长度。下面是带边界处理的写法：

coefs, freqs = pywt.cwt( x, scales, 'cmor1.5-1.0', sampling_period=1/fs, method='fft', # FFT 方式计算，速度更快 mode='periodization' ) # 丢弃首尾各 2% 的数据点，避免边界宽带影响 trim = int(len(t) * 0.02) coefs = coefs[:, trim:-trim]

method='fft'适合对整段信号做连续分析，因为 FFT 假定信号周期延拓，所以要求mode也匹配周期延拓。如果信号不是严格周期，改用method='conv'会稍微慢一些，但边界响应更接近直接积分。参数调整后要用下一节的自测信号确认频率峰位置，不能让边界亮带影响判断。

5. 用“已知频率信号”校验 wavelet.rar 里的 CWT 参数

不管从 wavelet.rar 里拿到的是.mat、.csv还是.npy数据，落地的第一步应该是构造一个已知频率的仿真信号，跑同一套 CWT 参数，检查频谱峰值频率是否准确。这个自检能一次性暴露采样率设置、尺度范围、小波中心频率三类错误。

def validate_cwt(scales, wavelet, fs=1000.0, f0=50.0): t = np.arange(0, 1, 1/fs) x = np.sin(2 * np.pi * f0 * t) coefs, freqs = pywt.cwt(x, scales, wavelet, sampling_period=1/fs) ridge_idx = np.argmax(np.abs(coefs), axis=0) ridge_freq = freqs[ridge_idx] detected = np.median(ridge_freq) rel_error = abs(detected - f0) / f0 print(f"detected={detected:.2f} Hz, error={rel_error:.2%}") return coefs, freqs

如果rel_error超过 5%，先检查scales是否覆盖了f0；再确认wavelet字符串里有没有把带宽和中心频率写反。实际中常见的错误是cmor1.5-1.0写成cmor1-1.5，结果中心频率变成 1.5，所有频率读数等比偏移 50%。请对比公式验证fc / (scale * Ts)中使用的 fc。

验证通过后，再把同一组参数应用到wavelet.rar提供的真实数据。如果真实数据谱图在预期频段之外还有持续亮带，先用scipy.signal.detrend去掉线性趋势，再看窗函数mode是否用了periodization。

如果真实数据和仿真信号的采样率不同，不要直接复用 scales，先按 4.2 的公式重算。pywt.cwt不会因为你传了sampling_period就自动把 scales 缩放到目标频段，它只会原样使用你传入的尺度数组。所以每次换数据先执行scale_min = fc * fs / fmax和scale_max = fc * fs / fmin这两行代码，再生成新的np.geomspace。检查完频率轴，再去看相位或包络，就不会被彩色图上的边界噪声干扰了。

本文还有配套的精品资源，点击获取