简介在时频分析中Morlet 小波由高斯窗与复指数构成能够兼顾时间与频率分辨率的平衡常用于生物医学、地震数据和经济时间序列等非平稳信号的分析。这份代码包面向 MATLAB 与 Python 开发者提供了创建、定义和使用 Morlet 小波的完整示例覆盖连续小波变换从参数设置、小波生成到结果可视化的主要流程适合信号处理初学者与需要快速实现时频分析的科研人员。包内共 6 个文件包含 2 个 MATLAB 脚本、1 个 Jupyter Notebook、1 个 PDF 原理说明、1 个 Markdown 文档和 1 个数据文件压缩包仅 3.21MB结构紧凑清晰便于直接阅读和运行。配套数据与说明文档有助于理解中心频率和尺度参数的选择逻辑也能帮助读者将方法迁移到自己的项目中。目前已有 634 人学习下载是一份轻量实用的时频分析入门参考。 前阵子处理一组轴承振动数据时我遇到一个典型问题信号里既有持续的周期振动又有短暂的非平稳冲击短时傅里叶变换的窗长怎么选都顾此失彼。后来换成复数 Morlet 小波做时频分析两个分量一下子分开了。当时先用 MATLAB 做了原型验证后面团队要跨平台复现又迁到 Python 里重写结果发现两边的定义、标注和归一化差异比想象中大得多。这篇内容就把 Morlet 小波在时域和频域中的创建、定义和使用讲清楚MATLAB 和 Python 都给出可直接运行的代码顺便把这两年在两个平台之间来回迁移时踩过的坑一并整理出来。1. 为什么偏偏是Morlet傅里叶处理不了的信号它能搞定1.1 傅里叶变换的“全局视野”局限傅里叶变换有一个很致命的问题它把整段信号拆成一系列正弦波的叠加得到的是“整个时间段里有哪些频率”却完全无法回答“某个频率是在什么时刻出现的”。对于平稳信号这完全够用但对于振动冲击、语音突变、脑电瞬态这类非平稳信号频谱图上只能看到一个平均结果时间信息彻底丢失。短时傅里叶变换STFT试图补救这个问题做法是把信号切成很多小段再对每一段做傅里叶变换。但一旦动手你就发现窗口的宽度是个两难选择窗口越长频率分辨率越高但时间上越来越“糊”窗口越短能定位瞬态位置但频率轴上相邻两个分量很难分开。我在实际分析中经常被这个矛盾卡住尤其是信号同时包含 30Hz 低频周期成分和 80Hz 高频短脉冲时想在一张图上同时看清楚两者几乎不可能。1.2 复值小波的独特位置小波变换的思路完全不同它不固定窗口长度而是用一系列“可伸缩”的小波函数去匹配信号。Morlet 小波属于其中非常特殊的一类它本质是一个复正弦波乘以高斯窗。因为带复数性质它能同时给出时频投影的幅值和相位这是分析瞬时频率、相位突变、振动模态等场景特别需要的性质。这么说吧如果你关心的是“能量在何时何频出现”用小波变换的模值就够了如果你还想进一步看“某个频率成分的相位如何随时间变化”那必须用复数小波Morlet 在工程和科研里几乎是默认选择。实际项目中我处理旋转机械的启停过程信号、语音信号共振峰追踪、脑电的时频能量分布时Morlet 都表现得很稳代码在 MATLAB 和 Python 里也都有很成熟的落地方案。1.3 什么时候该换其他小波Morlet 不是万能的。它频率选择性好对平滑的周期分量非常友好但对尖锐瞬态的定位不如 Mexican HatDOG 小波或 Haar 小波。如果信号里主要是大幅值的短暂冲击Morlet 画出来的时频图会在冲击位置拉出一条很宽的频带视觉上容易误判。我的经验是先看任务目标——要追相位、做瞬时频率估计优先 Morlet要检测突变点、边缘定位换一阶或二阶 DOG 小波更合适。这个选择在 MATLAB 和 Python 里都是一样的逻辑。2. Morlet小波的数学骨架写代码前必须看懂的三个参数2.1 时域定义复数 Morlet 小波的常用定义长这样ψ(t) (1 / √(π·f_b)) · exp(2jπ·f_c·t) · exp(-t² / f_b)拆开看它就是三个部分的乘积第一项是归一化系数保证小波能量基本恒定第二项是中心频率为 f_c 的复正弦负责“选频”第三项是标准差受 f_b 控制的高斯窗负责“定位时间”。这里两个参数最需要理解清楚f_c 是小波中心频率决定小波在频域里瞄准哪个频段f_b 是带宽参数直接决定高斯窗的展宽f_b 越大时域波形越宽频域带宽越窄频率分辨率越高但时间定位能力下降。举个直观例子我调参数时发现f_b 从 1 改成 4 之后时频图上的频率脊明显变细但时间方向上瞬态起止位置变得模糊。这个现象就是海森堡不确定性原理在实际处理中的体现你在代码里改参数本质上就在时间分辨率和频率分辨率之间做权衡。2.2 频域形式与尺度-频率换算Morlet 小波的傅里叶变换依然是一个高斯函数中心位于 f_c 附近。正是因为它频域上是个“窄带”高斯包络所以能很好地挑选特定频率。实际使用时我们很少直接用连续时间 t 去卷而是通过“尺度”去控制小波伸缩f_actual f_c / scale × fs这里的 fs 是采样率scale 是连续小波变换里的尺度参数。比如 fs 1000Hzf_c 1.0scale 20那实际对应频率就是 1000 / 20 50Hz。注意这个换算在 MATLAB 内置函数和 Python 工具包里都成立只是 API 的传参方式不同后面我会分别演示。2.3 归一化方式不同结果数值差很多写代码前还要做好一个心理准备MATLAB 和 Python 各工具包对 Morlet 小波采用的归一化方式不完全一致即使同样输入一段信号得到的小波系数幅值也可能差一个比例因子。这不代表谁算错了只代表它们对小波函数能量归一化的约定不同。在做跨平台结果对照时我的做法是先输入一个已知幅值的正弦信号分别跑一遍两边代码用峰值幅度或能量比值求出校正系数再统一量纲。直接拿两边的 raw 系数比大小没有任何意义这点后面实战部分会再演示。3. MATLAB实现从手工构造到cwt内置函数3.1 手工构造Morlet并验证频域形状在 MATLAB 里手工构造 Morlet 小波不复杂关键是别把时间轴算错。假设采样率 fs 1000中心频率 fc 10带宽 fb 1时间范围取 -2 到 2 秒fs 1000; fc 10; fb 1; t -2:1/fs:2; % 复数 Morlet omega 2 * pi * fc; psi_t 1/sqrt(pi*fb) * exp(1i*omega*t) .* exp(-(t.^2)/fb); % 频域验证 f linspace(-fs/2, fs/2, length(t)); psi_f fftshift(fft(psi_t)); plot(f, abs(psi_f)); xlabel(Frequency (Hz)); ylabel(|Psi(f)|);跑完之后你会看到频域幅值谱在 10Hz 附近形成一个漂亮的高斯峰。这里有个细节如果不用 1i 而是用 i 当虚数单位一旦前面重新赋值过 i 就很容易出隐晦的 bug所以我习惯全程使用 1i。这个小习惯在长时间脚本里能省不少排查时间。3.2 用内置cwt函数做时频分析R2016b 之后MATLAB 的 cwt 函数接口比老版本友好很多直接传原始信号和采样率即可% 生成测试信号 fs 1000; t 0:1/fs:1; x 0.6 * sin(2*pi*30*t) .* (t 0.2 t 0.5) ... 0.8 * sin(2*pi*80*t) .* (t 0.6 t 0.9); % 连续小波变换默认使用解析 Morlet [wt, freq] cwt(x, fs); % 绘制时频图 figure; imagesc(t, freq, abs(wt)); axis xy; xlabel(Time (s)); ylabel(Frequency (Hz));新版本 cwt 返回的 freq 已经是实际频率轴画图直接对应到 Hz非常方便。它内部默认用的就是解析 Morlet 小波所以做通用时频分析时通常不需要自己手工构造小波再去做卷积。3.3 MATLAB的版本坑与参数对齐老版本R2016a 及之前的 cwt 是另一套语法cwt(x, scales, morl)频率轴要用 scal2frq(scales, morl, fs) 或 centfrq 手动换算。如果你在网上搜到一堆旧代码一对比会发现新老版本返回结果的行列方向都不同new cwt 返回的是 time × frequency 的矩阵频率行、时间列而老版本是 scale × time 的矩阵。这个差异很容易让人在画图时想当然地转置结果时频图看上去完全错乱。我的建议是新项目直接用新接口维护老代码时先确认 cwt 的返回维度再决定是否 transpose。另外新接口的 cwt 默认自带边界锥形区域cone of influence信息画图时能直观看到哪些区域受边界效应影响我一般只读有效范围内的系数。4. Python实现scipy和PyWavelets两条路线4.1 用scipy手工构造与signal.cwtPython 里最基础的路线是复用 numpy 构造 Morlet再用 scipy.signal.cwt 做变换。手工构造小波时代码和 MATLAB 几乎一样import numpy as np from scipy import signal import matplotlib.pyplot as plt fs 1000 fc 10 fb 1 t np.arange(-2, 2, 1/fs) omega 2 * np.pi * fc psi_t 1/np.sqrt(np.pi*fb) * np.exp(1j*omega*t) * np.exp(-(t**2)/fb) # 频域验证 f np.fft.fftshift(np.fft.fftfreq(len(t), 1/fs)) psi_f np.fft.fftshift(np.fft.fft(psi_t)) plt.plot(f, np.abs(psi_f))scipy.signal.cwt 需要你传一个母小波函数和一组宽度 scales。它内部会按 scale 对母小波进行伸缩再和信号做卷积。一个常见写法是scales np.arange(1, 100) coefs signal.cwt(x, signal.morlet2, scales, s1)这里有一个必须注意的点signal.morlet2 的 s 参数不是中心频率而是与尺度配合的窗口宽度参数。它的定义是 exp(2jπx)·exp(-0.5·(x/s)²)默认不带能量归一化所以算出来的系数量级和 MATLAB cwt 结果完全不同不要拿数值直接对比。4.2 用pywt.cwt一步到位更省事的方案是用 PyWaveletsimport pywt scales np.arange(1, 100) coefs, freqs pywt.cwt(x, scales, cmor1.0-1.0, sampling_period1/fs) plt.pcolormesh(t, freqs, np.abs(coefs), shadinggouraud) plt.ylabel(Frequency (Hz)) plt.xlabel(Time (s))这里的 cmor1.0-1.0 是复数 Morlet 的小波族命名格式前一个数是带宽参数 f_b后一个数是中心频率 f_c。比如 cmor1.5-1.0 就是 fb1.5、fc1.0。pywt.cwt 支持直接传 sampling_period返回的 freqs 就是实际频率这个设计比 scipy 那条路顺手很多。实际使用时scales 的选择直接决定频率轴覆盖范围。经验公式是最小尺度对应最高分析频率最大尺度对应最低分析频率。如果你用 fc1.0 的 cmor且采样率 fs1000那么 frequency fs / scale。想看 5Hz 到 100Hz 的频段scale 取 10 到 200 即可。4.3 和MATLAB结果对应的关键点跨平台复现时我最常被问到的是“为什么两边画出来的图趋势一样数值差别很大”。原因前面说过归一化不同。pywt 的 cwt 默认会做 L1 归一化而 MATLAB 的 cwt 内部又有自己的一套约定两边系数直接差的不是一个固定倍数还可能与尺度相关。我在处理项目报告时会把两边的时频图都改成相对幅值比如把整张图的最大值归一化为 1这样能耗分布、瞬时频率脊的位置完全一致就不用来回解释数值差异了。如果你非要让两边绝对幅值对齐那就准备一个已知幅值正弦信号做校准算出逐尺度的比例系数再套用到目标信号上。这个方法不优雅但实测有效。5. 实测非平稳信号分离与调试避坑总结5.1 构造一个多分量验证信号纸上谈兵没意思我用一个典型的多分量信号来走完整流程。信号包含三个成分50Hz 持续正弦、120Hz 的短时冲击0.3 到 0.5 秒以及线性调频成分从 60Hz 扫到 90Hzfs 1000 t np.arange(0, 1, 1/fs) x (0.8 * np.sin(2*np.pi*50*t) 1.2 * np.sin(2*np.pi*120*t) * ((t 0.3) (t 0.5)) 0.5 * np.sin(2*np.pi*(60 30*t)*t))这个信号用普通傅里叶变换只能看到三个模糊的峰完全分不清 120Hz 的冲击到底发生在什么时候。用 pywt 的 cmor1.5-1.0 小波做变换后三个分量的时频轨迹非常清楚50Hz 是一条持续的直线120Hz 是时间窗口内的一条短线60-90Hz 是斜向上的带状轨迹三者相互独立能量互不串扰。这段测试是我最常用的验证方式用来检查小波参数和尺度范围是否设置正确。任何新环境配置好后先跑这段信号如果时频图能干净分出三个分量说明基础链路没问题如果糊成一团十有八九是频率轴换算错了。5.2 边界锥形区与边缘假影看时频图时最容易误导人的是图像左右两边的异常亮带。这不是信号真的在那段时间有强能量而是小波在信号边界处只覆盖了部分窗口导致卷积结果边缘被截断产生假影。pywt 和 MATLAB 内部虽然都会做边界延拓但无法完全消除这个问题。判断有效区域有个标准做法从左右两端各去掉小波最大尺度对应的时间跨度。实际操作中我一般直接用 MATLAB cwt 自带的锥形影响区COI显示范围或者手动在 Python 里按最大尺度的 1/2 时间长度裁剪图像。尺度范围越大边界效应影响越宽这个代价往往容易被新手忽略。5.3 三个花时间最多的坑第一个坑是尺度范围设置不合理。有人直接拿 1 到 10000 去跑 pywt结果频率轴几乎全挤在低频区图像完全没法看。正确做法是先敲定你关心的频段再用 fs / f 反推 scale 范围。第二个坑是忘了检查小波的中心频率。cmor 的 f_c 不同同样的 scale 对应的实际频率就不同。比如 cmor1.5-0.5 的中心频率是 0.5那么频率换算要按实际中心频率计算不能想当然套 fc1 的公式。第三个坑是时频图纵轴方向搞反。MATLAB 默认 imagesc 的 Y 轴是自上而下递增的如果不加 axis xy 或者不配合 pcolormesh画出来的图频率范围是反的。这种低级错误不看坐标轴刻度很难发现我见过不止一个同事对着倒置的时频图分析半天最后发现只是显示问题。6. 两个平台迁移时可以直接照抄的经验6.1 代码结构建议我在工程里把时频分析封装成几个独立模块小波构造、时频变换、绘图、边界掩码。MATLAB 和 Python 各写一套但接口保持一致输入都是 signal、fs、f_c、f_b、频率范围输出都是 time、freq、coefs这样切换平台时只有底层实现变化上层逻辑完全不用动。这个小设计看着多余但在项目从原型到产品迭代时省了非常多时间。6.2 常用参数参考表目标推荐 f_c推荐 f_b备注通用时频分析1.01.0~2.0平衡时间与频率分辨率瞬时频率估计1.04~8频率分辨率优先瞬态冲击检测1.00.5~1.0时间定位优先频带会变宽低频缓变振动0.52~4避免尺度范围过大这个表不是标准答案只是我从多次屡试不爽的调试里沉淀出的起点参数。实际信号类型变了参数应该以时频图效果为准再做微调。6.3 调参顺序和判断方法我调整参数有一套固定顺序先定频率分析范围再根据范围反推尺度然后调 f_b 控制频率脊宽度最后裁剪边界。每改一次参数就重新看一遍测试信号的时频图用已知成分来验证效果。如果人工信号都分不清就不要指望实际信号能有什么奇迹。还有个小技巧小波系数幅值不能直接理解成原信号幅值它是经过小波函数卷积后的投影强度。需要提取具体分量波形时不要试图从小波系数里直接还原原始信号正确做法是先用时频图确定分量所在的时频区间再用滤波器或傅里叶逆变换重建该频段的时域波形。这个弯路我走过一次耗时两天最后才明白小波系数量值和原信号幅度之间不是简单映射关系。Morlet 小波在 MATLAB 和 Python 里实现都不难难的是理解每个参数背后的物理含义和两套函数库的约定差异。只要把数学定义、尺度频率换算和边界效应这三件事先想清楚剩下的基本都是体力活。本文还有配套的精品资源点击获取 Morlet小波时频分析实战：MATLAB与Python实现全解析
简介在时频分析中Morlet 小波由高斯窗与复指数构成能够兼顾时间与频率分辨率的平衡常用于生物医学、地震数据和经济时间序列等非平稳信号的分析。这份代码包面向 MATLAB 与 Python 开发者提供了创建、定义和使用 Morlet 小波的完整示例覆盖连续小波变换从参数设置、小波生成到结果可视化的主要流程适合信号处理初学者与需要快速实现时频分析的科研人员。包内共 6 个文件包含 2 个 MATLAB 脚本、1 个 Jupyter Notebook、1 个 PDF 原理说明、1 个 Markdown 文档和 1 个数据文件压缩包仅 3.21MB结构紧凑清晰便于直接阅读和运行。配套数据与说明文档有助于理解中心频率和尺度参数的选择逻辑也能帮助读者将方法迁移到自己的项目中。目前已有 634 人学习下载是一份轻量实用的时频分析入门参考。 前阵子处理一组轴承振动数据时我遇到一个典型问题信号里既有持续的周期振动又有短暂的非平稳冲击短时傅里叶变换的窗长怎么选都顾此失彼。后来换成复数 Morlet 小波做时频分析两个分量一下子分开了。当时先用 MATLAB 做了原型验证后面团队要跨平台复现又迁到 Python 里重写结果发现两边的定义、标注和归一化差异比想象中大得多。这篇内容就把 Morlet 小波在时域和频域中的创建、定义和使用讲清楚MATLAB 和 Python 都给出可直接运行的代码顺便把这两年在两个平台之间来回迁移时踩过的坑一并整理出来。1. 为什么偏偏是Morlet傅里叶处理不了的信号它能搞定1.1 傅里叶变换的“全局视野”局限傅里叶变换有一个很致命的问题它把整段信号拆成一系列正弦波的叠加得到的是“整个时间段里有哪些频率”却完全无法回答“某个频率是在什么时刻出现的”。对于平稳信号这完全够用但对于振动冲击、语音突变、脑电瞬态这类非平稳信号频谱图上只能看到一个平均结果时间信息彻底丢失。短时傅里叶变换STFT试图补救这个问题做法是把信号切成很多小段再对每一段做傅里叶变换。但一旦动手你就发现窗口的宽度是个两难选择窗口越长频率分辨率越高但时间上越来越“糊”窗口越短能定位瞬态位置但频率轴上相邻两个分量很难分开。我在实际分析中经常被这个矛盾卡住尤其是信号同时包含 30Hz 低频周期成分和 80Hz 高频短脉冲时想在一张图上同时看清楚两者几乎不可能。1.2 复值小波的独特位置小波变换的思路完全不同它不固定窗口长度而是用一系列“可伸缩”的小波函数去匹配信号。Morlet 小波属于其中非常特殊的一类它本质是一个复正弦波乘以高斯窗。因为带复数性质它能同时给出时频投影的幅值和相位这是分析瞬时频率、相位突变、振动模态等场景特别需要的性质。这么说吧如果你关心的是“能量在何时何频出现”用小波变换的模值就够了如果你还想进一步看“某个频率成分的相位如何随时间变化”那必须用复数小波Morlet 在工程和科研里几乎是默认选择。实际项目中我处理旋转机械的启停过程信号、语音信号共振峰追踪、脑电的时频能量分布时Morlet 都表现得很稳代码在 MATLAB 和 Python 里也都有很成熟的落地方案。1.3 什么时候该换其他小波Morlet 不是万能的。它频率选择性好对平滑的周期分量非常友好但对尖锐瞬态的定位不如 Mexican HatDOG 小波或 Haar 小波。如果信号里主要是大幅值的短暂冲击Morlet 画出来的时频图会在冲击位置拉出一条很宽的频带视觉上容易误判。我的经验是先看任务目标——要追相位、做瞬时频率估计优先 Morlet要检测突变点、边缘定位换一阶或二阶 DOG 小波更合适。这个选择在 MATLAB 和 Python 里都是一样的逻辑。2. Morlet小波的数学骨架写代码前必须看懂的三个参数2.1 时域定义复数 Morlet 小波的常用定义长这样ψ(t) (1 / √(π·f_b)) · exp(2jπ·f_c·t) · exp(-t² / f_b)拆开看它就是三个部分的乘积第一项是归一化系数保证小波能量基本恒定第二项是中心频率为 f_c 的复正弦负责“选频”第三项是标准差受 f_b 控制的高斯窗负责“定位时间”。这里两个参数最需要理解清楚f_c 是小波中心频率决定小波在频域里瞄准哪个频段f_b 是带宽参数直接决定高斯窗的展宽f_b 越大时域波形越宽频域带宽越窄频率分辨率越高但时间定位能力下降。举个直观例子我调参数时发现f_b 从 1 改成 4 之后时频图上的频率脊明显变细但时间方向上瞬态起止位置变得模糊。这个现象就是海森堡不确定性原理在实际处理中的体现你在代码里改参数本质上就在时间分辨率和频率分辨率之间做权衡。2.2 频域形式与尺度-频率换算Morlet 小波的傅里叶变换依然是一个高斯函数中心位于 f_c 附近。正是因为它频域上是个“窄带”高斯包络所以能很好地挑选特定频率。实际使用时我们很少直接用连续时间 t 去卷而是通过“尺度”去控制小波伸缩f_actual f_c / scale × fs这里的 fs 是采样率scale 是连续小波变换里的尺度参数。比如 fs 1000Hzf_c 1.0scale 20那实际对应频率就是 1000 / 20 50Hz。注意这个换算在 MATLAB 内置函数和 Python 工具包里都成立只是 API 的传参方式不同后面我会分别演示。2.3 归一化方式不同结果数值差很多写代码前还要做好一个心理准备MATLAB 和 Python 各工具包对 Morlet 小波采用的归一化方式不完全一致即使同样输入一段信号得到的小波系数幅值也可能差一个比例因子。这不代表谁算错了只代表它们对小波函数能量归一化的约定不同。在做跨平台结果对照时我的做法是先输入一个已知幅值的正弦信号分别跑一遍两边代码用峰值幅度或能量比值求出校正系数再统一量纲。直接拿两边的 raw 系数比大小没有任何意义这点后面实战部分会再演示。3. MATLAB实现从手工构造到cwt内置函数3.1 手工构造Morlet并验证频域形状在 MATLAB 里手工构造 Morlet 小波不复杂关键是别把时间轴算错。假设采样率 fs 1000中心频率 fc 10带宽 fb 1时间范围取 -2 到 2 秒fs 1000; fc 10; fb 1; t -2:1/fs:2; % 复数 Morlet omega 2 * pi * fc; psi_t 1/sqrt(pi*fb) * exp(1i*omega*t) .* exp(-(t.^2)/fb); % 频域验证 f linspace(-fs/2, fs/2, length(t)); psi_f fftshift(fft(psi_t)); plot(f, abs(psi_f)); xlabel(Frequency (Hz)); ylabel(|Psi(f)|);跑完之后你会看到频域幅值谱在 10Hz 附近形成一个漂亮的高斯峰。这里有个细节如果不用 1i 而是用 i 当虚数单位一旦前面重新赋值过 i 就很容易出隐晦的 bug所以我习惯全程使用 1i。这个小习惯在长时间脚本里能省不少排查时间。3.2 用内置cwt函数做时频分析R2016b 之后MATLAB 的 cwt 函数接口比老版本友好很多直接传原始信号和采样率即可% 生成测试信号 fs 1000; t 0:1/fs:1; x 0.6 * sin(2*pi*30*t) .* (t 0.2 t 0.5) ... 0.8 * sin(2*pi*80*t) .* (t 0.6 t 0.9); % 连续小波变换默认使用解析 Morlet [wt, freq] cwt(x, fs); % 绘制时频图 figure; imagesc(t, freq, abs(wt)); axis xy; xlabel(Time (s)); ylabel(Frequency (Hz));新版本 cwt 返回的 freq 已经是实际频率轴画图直接对应到 Hz非常方便。它内部默认用的就是解析 Morlet 小波所以做通用时频分析时通常不需要自己手工构造小波再去做卷积。3.3 MATLAB的版本坑与参数对齐老版本R2016a 及之前的 cwt 是另一套语法cwt(x, scales, morl)频率轴要用 scal2frq(scales, morl, fs) 或 centfrq 手动换算。如果你在网上搜到一堆旧代码一对比会发现新老版本返回结果的行列方向都不同new cwt 返回的是 time × frequency 的矩阵频率行、时间列而老版本是 scale × time 的矩阵。这个差异很容易让人在画图时想当然地转置结果时频图看上去完全错乱。我的建议是新项目直接用新接口维护老代码时先确认 cwt 的返回维度再决定是否 transpose。另外新接口的 cwt 默认自带边界锥形区域cone of influence信息画图时能直观看到哪些区域受边界效应影响我一般只读有效范围内的系数。4. Python实现scipy和PyWavelets两条路线4.1 用scipy手工构造与signal.cwtPython 里最基础的路线是复用 numpy 构造 Morlet再用 scipy.signal.cwt 做变换。手工构造小波时代码和 MATLAB 几乎一样import numpy as np from scipy import signal import matplotlib.pyplot as plt fs 1000 fc 10 fb 1 t np.arange(-2, 2, 1/fs) omega 2 * np.pi * fc psi_t 1/np.sqrt(np.pi*fb) * np.exp(1j*omega*t) * np.exp(-(t**2)/fb) # 频域验证 f np.fft.fftshift(np.fft.fftfreq(len(t), 1/fs)) psi_f np.fft.fftshift(np.fft.fft(psi_t)) plt.plot(f, np.abs(psi_f))scipy.signal.cwt 需要你传一个母小波函数和一组宽度 scales。它内部会按 scale 对母小波进行伸缩再和信号做卷积。一个常见写法是scales np.arange(1, 100) coefs signal.cwt(x, signal.morlet2, scales, s1)这里有一个必须注意的点signal.morlet2 的 s 参数不是中心频率而是与尺度配合的窗口宽度参数。它的定义是 exp(2jπx)·exp(-0.5·(x/s)²)默认不带能量归一化所以算出来的系数量级和 MATLAB cwt 结果完全不同不要拿数值直接对比。4.2 用pywt.cwt一步到位更省事的方案是用 PyWaveletsimport pywt scales np.arange(1, 100) coefs, freqs pywt.cwt(x, scales, cmor1.0-1.0, sampling_period1/fs) plt.pcolormesh(t, freqs, np.abs(coefs), shadinggouraud) plt.ylabel(Frequency (Hz)) plt.xlabel(Time (s))这里的 cmor1.0-1.0 是复数 Morlet 的小波族命名格式前一个数是带宽参数 f_b后一个数是中心频率 f_c。比如 cmor1.5-1.0 就是 fb1.5、fc1.0。pywt.cwt 支持直接传 sampling_period返回的 freqs 就是实际频率这个设计比 scipy 那条路顺手很多。实际使用时scales 的选择直接决定频率轴覆盖范围。经验公式是最小尺度对应最高分析频率最大尺度对应最低分析频率。如果你用 fc1.0 的 cmor且采样率 fs1000那么 frequency fs / scale。想看 5Hz 到 100Hz 的频段scale 取 10 到 200 即可。4.3 和MATLAB结果对应的关键点跨平台复现时我最常被问到的是“为什么两边画出来的图趋势一样数值差别很大”。原因前面说过归一化不同。pywt 的 cwt 默认会做 L1 归一化而 MATLAB 的 cwt 内部又有自己的一套约定两边系数直接差的不是一个固定倍数还可能与尺度相关。我在处理项目报告时会把两边的时频图都改成相对幅值比如把整张图的最大值归一化为 1这样能耗分布、瞬时频率脊的位置完全一致就不用来回解释数值差异了。如果你非要让两边绝对幅值对齐那就准备一个已知幅值正弦信号做校准算出逐尺度的比例系数再套用到目标信号上。这个方法不优雅但实测有效。5. 实测非平稳信号分离与调试避坑总结5.1 构造一个多分量验证信号纸上谈兵没意思我用一个典型的多分量信号来走完整流程。信号包含三个成分50Hz 持续正弦、120Hz 的短时冲击0.3 到 0.5 秒以及线性调频成分从 60Hz 扫到 90Hzfs 1000 t np.arange(0, 1, 1/fs) x (0.8 * np.sin(2*np.pi*50*t) 1.2 * np.sin(2*np.pi*120*t) * ((t 0.3) (t 0.5)) 0.5 * np.sin(2*np.pi*(60 30*t)*t))这个信号用普通傅里叶变换只能看到三个模糊的峰完全分不清 120Hz 的冲击到底发生在什么时候。用 pywt 的 cmor1.5-1.0 小波做变换后三个分量的时频轨迹非常清楚50Hz 是一条持续的直线120Hz 是时间窗口内的一条短线60-90Hz 是斜向上的带状轨迹三者相互独立能量互不串扰。这段测试是我最常用的验证方式用来检查小波参数和尺度范围是否设置正确。任何新环境配置好后先跑这段信号如果时频图能干净分出三个分量说明基础链路没问题如果糊成一团十有八九是频率轴换算错了。5.2 边界锥形区与边缘假影看时频图时最容易误导人的是图像左右两边的异常亮带。这不是信号真的在那段时间有强能量而是小波在信号边界处只覆盖了部分窗口导致卷积结果边缘被截断产生假影。pywt 和 MATLAB 内部虽然都会做边界延拓但无法完全消除这个问题。判断有效区域有个标准做法从左右两端各去掉小波最大尺度对应的时间跨度。实际操作中我一般直接用 MATLAB cwt 自带的锥形影响区COI显示范围或者手动在 Python 里按最大尺度的 1/2 时间长度裁剪图像。尺度范围越大边界效应影响越宽这个代价往往容易被新手忽略。5.3 三个花时间最多的坑第一个坑是尺度范围设置不合理。有人直接拿 1 到 10000 去跑 pywt结果频率轴几乎全挤在低频区图像完全没法看。正确做法是先敲定你关心的频段再用 fs / f 反推 scale 范围。第二个坑是忘了检查小波的中心频率。cmor 的 f_c 不同同样的 scale 对应的实际频率就不同。比如 cmor1.5-0.5 的中心频率是 0.5那么频率换算要按实际中心频率计算不能想当然套 fc1 的公式。第三个坑是时频图纵轴方向搞反。MATLAB 默认 imagesc 的 Y 轴是自上而下递增的如果不加 axis xy 或者不配合 pcolormesh画出来的图频率范围是反的。这种低级错误不看坐标轴刻度很难发现我见过不止一个同事对着倒置的时频图分析半天最后发现只是显示问题。6. 两个平台迁移时可以直接照抄的经验6.1 代码结构建议我在工程里把时频分析封装成几个独立模块小波构造、时频变换、绘图、边界掩码。MATLAB 和 Python 各写一套但接口保持一致输入都是 signal、fs、f_c、f_b、频率范围输出都是 time、freq、coefs这样切换平台时只有底层实现变化上层逻辑完全不用动。这个小设计看着多余但在项目从原型到产品迭代时省了非常多时间。6.2 常用参数参考表目标推荐 f_c推荐 f_b备注通用时频分析1.01.0~2.0平衡时间与频率分辨率瞬时频率估计1.04~8频率分辨率优先瞬态冲击检测1.00.5~1.0时间定位优先频带会变宽低频缓变振动0.52~4避免尺度范围过大这个表不是标准答案只是我从多次屡试不爽的调试里沉淀出的起点参数。实际信号类型变了参数应该以时频图效果为准再做微调。6.3 调参顺序和判断方法我调整参数有一套固定顺序先定频率分析范围再根据范围反推尺度然后调 f_b 控制频率脊宽度最后裁剪边界。每改一次参数就重新看一遍测试信号的时频图用已知成分来验证效果。如果人工信号都分不清就不要指望实际信号能有什么奇迹。还有个小技巧小波系数幅值不能直接理解成原信号幅值它是经过小波函数卷积后的投影强度。需要提取具体分量波形时不要试图从小波系数里直接还原原始信号正确做法是先用时频图确定分量所在的时频区间再用滤波器或傅里叶逆变换重建该频段的时域波形。这个弯路我走过一次耗时两天最后才明白小波系数量值和原信号幅度之间不是简单映射关系。Morlet 小波在 MATLAB 和 Python 里实现都不难难的是理解每个参数背后的物理含义和两套函数库的约定差异。只要把数学定义、尺度频率换算和边界效应这三件事先想清楚剩下的基本都是体力活。本文还有配套的精品资源点击获取