简介面向小波模极大值方法在MATLAB中的实现需求这份资源以单脚本形式提供了完整的算法演示适用于信号检测、特征提取与图像分析方向的学习者以及需要处理突变信号或噪声背景的研究人员。压缩包内仅含1个m文件整体体积仅563B但脚本流程清晰依次涵盖小波分解、各尺度模极大值计算与定位、阈值去噪处理以及基于极大值点的信号重构或特征输出可直接在MATLAB中运行修改。已有227人学习了该资源适合希望理解“cgau”连续高斯小波特性、掌握wavedec与waverec函数配合使用的读者尤其对想快速入门小波模极大值原理的初学者颇具参考价值。通过研读并运行该脚本读者既能观察到信号突变点的定位效果也能借鉴其代码框架将小波模极大值方法迁移至地震数据、医学成像或金融时间序列等实际场景中是一份轻量而实用的动手学习资料。1. 小波模极大值检测奇异点先看信号在哪一刻变了把压缩包名称拆开“qiyidian”是“奇异点”的拼音“小波模极大值”才是核心算法信号在某一瞬间发生跳变、脉冲或斜率突变时普通差分和FFT只能告诉你“有大变化”却回答不了变化发生在第几个采样点、突变有多尖锐。小波模极大值方法把原始信号放到多尺度空间里观察每个突变点的小波系数模极大值会沿尺度方向形成一条极值线其衰减速率由Lipschitz指数决定。把这个指数估算出来就能把“突变”从定性变成定量。下面这套流程面向故障诊断、信号处理和MATLAB工程实践从算法原理、代码实现、参数调整到验证技巧完整走一遍从原始波形到奇异点列表的路线。新手能按步骤复现熟手也能看到噪声条件下参数取舍的边界。2. 小波模极大值的数学基础Lipschitz指数与极值线传播2.1 Lipschitz指数度量奇异强度奇异点不是非黑即白的“跳变”同样的电压突降电网上相邻节点捕捉到的陡峭程度不同同样的边缘像素经图像退化后梯度也会变缓。Lipschitz指数把“陡峭程度”收敛成一个实数。若x(t)在t0附近满足|x(t) - x(t0)| ≤ K·|t - t0|^α则称x(t)在t0处具有Lipschitz指数α。α1对应光滑可微点α0对应有界间断阶跃α-1对应比阶跃更尖锐的冲击理想脉冲0和1之间的值对应斜坡状突变。指数越小奇异越强。在离散信号里直接算这个指数很困难因为采样点落在突变处的相位会严重影响差分结果。小波变换的优势在于它相当于用可伸缩的窗口对信号做带通卷积每个尺度的输出都携带该尺度下突变轮廓的信息。Mallat证明了小波系数模极大值满足log|Wx(s,t)| ≤ log K α·log s这给了工程上一个可操作的估计办法沿极值线取多个尺度的小波系数模在双对数坐标里拟合直线斜率即α。实现这一步所需的全部工作就是先把小波变换系数算对再把极值点从二维矩阵里准确拎出来。2.2 小波基的选择影响极值线含义先将数学结论落到小波基选择上。用于奇异点检测的小波必须满足两个条件具有一阶消失矩或以上能够反映信号的局部变化尽可能平滑避免小波自身的高频抖动在模极大值提取阶段产生伪点。最常用的是高斯函数的一阶导数即gauss1小波。gauss1在时域是单峰形态的疏波奇数对称模板本身没有直流分量卷积结果近似于对信号做一阶微分后再平滑。不同小波基对应不同的“极值位置含义”。表2-1列出的定位关系决定后续追踪逻辑建议在写代码前先想清楚自己的目标是找信号本身突变用gauss1还是找信号拐点用gauss2。表2-1 不同小波基对奇异点检测的定位差异小波基模极大值对应定位精度典型用途haar小波差分突变点一般快速粗检测gauss1信号一阶导数过零点高奇异点检测首选gauss2信号拐点中图像边缘检测morse时频脊线点低时频分析不推荐做奇异性估计为什么不建议直接用MATLAB的cwt默认morse小波估Lipschitz指数morse小波在尺度频带上的定位更偏向时频分析其相位特性和极值衰减规律比gauss1复杂拟合出的斜率波动大解释起来也更困难。gauss1的另一个好处是模板解析式简单便于在调试时手算特定尺度下的期望模值。如果需要确认当前MATLAB版本是否支持直接调用gauss1在命令行执行waveinfo(gaus)查看即可。2.3 极值线在尺度空间中的传播与偏移一个小窍门若在最小尺度上检测到模极大值并且随着尺度增大该极大值点持续存在、位置漂移不超过各自尺度下小波支撑半径则基本可判定这是一个真实奇异点。理论来源是奇异点周围小波系数构成锥形影响区极值线穿过该区域一直延伸到最大尺度。反过来噪声形成的极值线覆盖尺度很短通常在2到3个尺度就消失。定位和偏移的量级需要估计。尺度s下gauss1的等效支撑半径大约为2.5s个采样间隔若该尺度下有另一个同极性奇异点在附近两条极值线会在大尺度方向汇合。实际检测时时间坐标取极值线最细尺度端的位置而不是把整条线平均这样能避免大尺度偏移引入的定位误差。下面这段代码是极值线追踪和Lipschitz拟合的基础也解释了为什么会用双对数坐标scales 2 .^ (1:0.5:7); % 指数增长的尺度序列 amp 0.8 * scales .^ (-0.35); % 假设α-0.35 loglog(scales, amp, o-); xlabel(尺度 s); ylabel(模极大值 |W|);若模值随尺度呈幂律衰减在双对数坐标里就是一条直线直线斜率即α。写代码之前先造出这样一组数据做基准实现阶段会少走很多弯路。3. MATLAB实现从连续小波变换到模极大值提取3.1 构造测试信号与gauss1小波核先搭建一个最小可运行流程。测试信号应包含至少两种尺度的突变阶跃和冲击。阶跃突变对应α0冲击对应α≈-1二者在模极大值提取中的行为差异足够明显。Fs 1000; t (0:999) / Fs; x sin(2*pi*50*t); % 干净背景 x(300:end) x(300:end) 0.6; % 阶跃突变位于第300点 x(700) x(700) 1.0; % 冲击突变位于第700点信号构造本身很简单第300点开始整体抬高形成阶跃第700点是单点脉冲。测试信号的作用是给整个流程提供一个可对照的真值检测完成后能直接拿检测位置与构造位置做差定位误差一目了然。接下来生成gauss1小波核。常见做法是先固定时间轴分辨率dt再按尺度s决定核长度。这样生成的每个尺度核都覆盖了高斯函数的完整支撑范围不会因为核太短把高频细节切掉dt 1 / Fs; nmax ceil(8 * s / dt); % 尺度s下核的半长度 nu -nmax:nmax; psi -(nu * dt / s) .* exp(-(nu * dt / s).^2 / 2); W conv(x, psi, same) * (dt / sqrt(s));这段代码的核心是把连续小波定义离散化小波模板ψ(u)的自变量u取为(n·dt)/s尺度越大核越长对应观察的波形窗口越宽。末尾乘上dt/√s是把积分写成黎曼和后的小波能量归一化保证不同尺度下的模值处于同一量级后面Lipschitz拟合才有意义。完整循环就是逐尺度做卷积把每一层的输出作为一行堆成W矩阵。信号两端用conv的same模式会默认补零制造边界假极值。建议先对x做镜像延拓xext wextend(sym, 2, x, floor(nmax/2)2); Wext conv(xext, psi, same); W Wext(floor(nmax/2)3 : end - floor(nmax/2)-2);延拓后出现假极值的位置会被裁剪掉这是数字实现里最容易忽略的一层也是“同样代码别人结果好、自己结果乱”的最常见原因。完整代码中需要把延拓放在尺度循环外一次完成避免每个尺度重复延拓造成边界不一致。3.2 在每个尺度行上提取模极大值候选点小波系数矩阵建好后奇异点检测的下一步是把每个尺度上|W|的局部极大值抽出来。这里要区分两个概念单尺度内的局部极大值和跨尺度同一条极值线上的点。前者用一维滑动比较或islocalmax即可threshold 0.1 * max(abs(W(:))); cand cell(length(scales), 1); for k 1:length(scales) w abs(W(k, :)); lm islocalmax(w, MinProminence, 0.05 * max(w)); lm lm (w threshold); pos find(lm); cand{k} [pos(:), w(lm(:))]; % [位置模值] endislocalmax的MinProminence参数值得展开说。没有这个参数时只要某点两侧比邻点大就会判定为极大值含噪信号里会出现大量假峰。MinProminence设为当前尺度最大模值的5%能滤掉明显的浅峰。threshold阈值再兜底一次把整体模值很低的候选点全部去掉。如果你的MATLAB版本较老没有islocalmax可以用diff(w)完成等价的局部极大值判定核心逻辑不变。cand{k}里存放的是二维数组第一列是时间位置第二列是对应的模值。后续极值线追踪完全依赖这两个字段所以提取阶段的过滤宁可偏严也不要偏松否则追踪代码会因为这些太密的候选点而连线错乱。3.3 跨尺度极值线追踪追踪的目标是把同一奇异点在不同尺度产生的模极大值连成线。常见做法是从最小尺度出发向大尺度找数最近邻当前尺度的候选点pos在上一尺度所有候选点中找时间差最小的点若差值小于允许范围就判定为同一条线。允许范围一般取0.2倍当前尺度对应的小波支撑宽度或干脆固定为3到5个采样点dtmax 0.2 * scales(k) * dt 2 * dt; ids zeros(size(cand{k}, 1), 1); for k 2:length(scales) c1 cand{k-1}; c2 cand{k}; for i 1:size(c1, 1) [mind, j] min(abs(c2(:, 1) - c1(i, 1))); if mind dtmax id assign_line(c1(i,:), c2(j,:)); % 把两点并到一条线 end end end此处的assign_line是你要自己维护的线id分配函数在正式工程里可以用图论里的连通分量来做避免多条线在赋值时互相覆盖。更简单的实现是把lines存成元胞数组每条线保存posSeq和ampSeq匹配成功就在线尾追加节点如果当前候选点与任何已有线都不匹配则新开一条线。追踪方向也可以从大尺度向小尺度反向进行但对噪声信号从小尺度出发更容易保持定位精度。追踪完成后一条有效的极值线至少要覆盖60%以上尺寸。这个比例是过滤噪声伪线的基本标准具体取值在第4章结合参数展开。3.4 用双对数拟合估计Lipschitz指数每条线的posSeq对应的尺度系数已知ampSeq由cand第二列得到。接下来只需要取对数后做线性拟合function alpha est_lipschitz(scales, amps, kRange) p polyfit(log(scales(kRange)), log(amps(kRange)), 1); alpha p(1); end调用例子如果某条线在尺度索引5到12之间均有数据就令kRange5:12得到该点的α估计值。为什么要舍弃两端的尺度最小尺度的模值容易受噪声扰动最大尺度的模值又可能被邻近奇异点通过极值线合并而干扰中间段更接近“只受该奇异点控制”的理论区域拟合稳定性最高。实际操作时通常做两次过滤第一次按覆盖率过滤线第二次按拟合残差过滤。若拟合残差的均方根超过0.3说明该线位置附近可能有两个互相影响的奇异点或者小波变换数值精度有问题需要回到第3.1节的延拓部分检查。拟合得到α后按表3-1做物理含义映射。表3-1 估计α与奇异类型的对应参考估计α范围奇异类型典型对象-1.2 ~ -0.8冲击局部放电脉冲、撞击信号-0.2 ~ 0.2阶跃相位切换、机械断裂0.3 ~ 0.7斜坡切变趋势拐点、缓慢饱和0.8 以上近似光滑一般不作为奇异点处理4. 参数选择与噪声场景下的小波模极大值调节4.1 尺度范围与尺度步长怎么定尺度范围的选择同时影响计算量和检测精度。最小尺度决定定位精度尺度1对应的gauss1核只覆盖约8个采样点定位误差也在1-2个采样点之内。如果信号采样率是1000Hz最小尺度取1就够若信号本身带宽有限例如传感器前端内置低通滤波则需要把最小尺度抬高到2-4避免检测出大量高频伪突变。最大尺度随着噪声水平和信号长度的变化而变化。理论上尺度取到信号长度的1/4以上没有意义因为此时整个信号基本都在小波支撑内奇异点已无法在时域上区分。经验上取信号长度的1/16到1/8。表4-1给出典型的初始参数组合实际使用时按信号能量分布微调最大尺度即可。表4-1 不同采样率下的尺度范围建议Fs信号长度最小尺度最大尺度步长128Hz2s1160.51kHz1s1640.510kHz0.5s21280.550kHz0.1s4641步长决定极值线匹配成功率。指数增长序列2^(a:b:c)中b取0.5时相邻尺度之间小波核长度相差约1.4倍核形状变化平缓追踪稳定。b取1时计算量小但相邻尺度间的模值可能突变对噪声敏感。推荐先从b0.5开始若信号很短或者只需要做粗检测再换到b1。4.2 阈值与显著性过滤固定百分比threshold0.1*max(|W|)在背景噪声很小时工作正常一旦噪声水平升高真实微弱奇异点的模值可能远低于背景峰值固定百分比就直接漏检。更稳健的做法是从中位数估计噪声能量noiseLevel median(abs(W(:))) / 0.6745; thr 4 * noiseLevel;0.6745是高斯噪声标准差与中位数绝对偏差的换算常数4倍噪声标准差对应的虚警率已经很低。该方法不需要知道信号的真实信噪比对工程应用更友好。如果信号中存在明显的有色噪声中位数会被低频成分抬高推荐改用每一尺度行单独计算中位数得到每层不同的阈值。噪声本身也会在尺度平面里形成极值点但这些极值点的显著性随尺度增大迅速衰退。若某条极值线覆盖尺度数少于总数一半其拟合出的α可信度非常低建议直接丢弃。阈值和覆盖率两个参数要一起调节只动一个往往压不住噪声。4.3 伪极值线的判别覆盖率与线长追踪完成后按线覆盖的尺度数过滤valid false(length(lines), 1); for i 1:length(lines) valid(i) numel(lines(i).amp) 0.6 * length(scales); end lines lines(valid);0.6这个比例需要结合尺度范围做微调。当尺度步长b0.5、总尺度数是13时真实极值线通常能覆盖全部13层噪声线一般只覆盖2到4层。如果最小尺度或最大尺度设置过短真线的覆盖率也会下降因此判断标准以“至少覆盖总尺度数的60%”为准不要直接写死为某个固定数量。白噪声在尺度谱上还有一个明显特征极值线分布均匀且密度高。过滤时可以把尺度平面的二维密度也纳入考量利用histcounts2统计单位区域的线数密度超标的区域整块剔除。这样处理对突发性强噪声更有效因为突发噪声通常只集中在某几个尺度上。4.4 低信噪比下的预平滑替代方案当信噪比降到10dB以下时只靠提高阈值往往压制不住噪声极值线需要在对小波系数做极值提取之前先做时间方向的中值滤波。常见做法是对每一尺度行w做medfilt1(w, m)m取该尺度小波主瓣宽度的一半。中值滤波在保边缘的同时能显著降低噪声峰缺点是小尺度端的极值位置会偏移几个样本点定位精度会略微下降。另一种方案是提高阈值的同时放宽dtmax。噪声让极值线位置抖动如果追踪窗口过紧真实极值线会被切碎成若干短线。此时可以把dtmax从固定3个采样点调整为0.4*scales(k)dt3dt。放宽追踪窗口的代价是邻近两个奇异点的极值线可能被误并为一条解决方法是检查线内模值的平滑性若同一条线内模值出现非单调反复则拆线重连。5. 验证技巧用合成信号检验小波模极大值的整条链路5.1 构造三类标准奇异点写完代码后第一件事不是跑真实数据而是用理论α已知的三类信号验证。阶跃点的Lipschitz指数理论为0脉冲点理论接近-1斜坡突变的理论在0.5附近。构造这样一支信号后运行前述函数比对检测位置与估计α。t (0:999) / 1000; x zeros(1, 1000); x(400) 1; % 脉冲 x(600:end) x(600:end) 1; % 阶跃 x(800:end) x(800:end) (0:199) / 200; % 斜坡切变脉冲在单点注入能量阶跃从第600点开始持久改变电平斜坡则在第800点开始引入斜率。三者的理论α分别接近-1、0和0.5足够覆盖多数故障信号中的突变形态。构造完成后分别调用极值提取和追踪函数观察三条极值线是否都出现。5.2 与理论指标做定量对照验证时要看的三个指标定位误差、α估计误差、极值线覆盖率。定位误差应小于最小尺度下的核半径α误差在±0.15以内算通过覆盖率应接近100%。任何一个不合格都优先回溯到边界延拓与阈值设置而不是直接改拟合代码。把这三个指标打印出来作为回归基准之后每改动一次参数跑一遍合成信号对比结果能快速发现是参数还是实现出了问题。若定位误差持续偏大可以在极值线追踪后用一次一维抛物线插值对峰值位置做亚像素修正这是不增加计算量就能提升定位精度的常用手段。5.3 按qiyidian.rar的组织习惯拆分函数标题里的qiyidian.rar让人联想到一个现成的MATLAB包。不管解压后原始文件怎么命名按“数据读取-尺度变换-极值提取-线追踪-指数估计”拆成五个函数去做以后换数据、换小波基时只需改对应层。例如换用gauss2只需要替换尺度变换里的模板函数其他模块完全不动。推荐的调用顺序W wt_sgram(x, scales, Fs); cand wt_findmax(W, MinProminence, 0.05); lines wt_track(cand, scales, dtmax, 0.3); result wt_lipschitz(lines, scales, MinCoverage, 0.6);参数名在函数接口里显式暴露比在脚本里写散落的魔数更容易维护。把合成信号验证脚本放到最前面作为一个独立的testSmoke.m形成每次改代码后必跑一次的回归习惯。这样即使后面接入的工程信号形态复杂也能保证小波模极大值这层核心算法始终处于可信状态。本文还有配套的精品资源点击获取
基于小波模极大值的信号奇异点检测与MATLAB实现
简介面向小波模极大值方法在MATLAB中的实现需求这份资源以单脚本形式提供了完整的算法演示适用于信号检测、特征提取与图像分析方向的学习者以及需要处理突变信号或噪声背景的研究人员。压缩包内仅含1个m文件整体体积仅563B但脚本流程清晰依次涵盖小波分解、各尺度模极大值计算与定位、阈值去噪处理以及基于极大值点的信号重构或特征输出可直接在MATLAB中运行修改。已有227人学习了该资源适合希望理解“cgau”连续高斯小波特性、掌握wavedec与waverec函数配合使用的读者尤其对想快速入门小波模极大值原理的初学者颇具参考价值。通过研读并运行该脚本读者既能观察到信号突变点的定位效果也能借鉴其代码框架将小波模极大值方法迁移至地震数据、医学成像或金融时间序列等实际场景中是一份轻量而实用的动手学习资料。1. 小波模极大值检测奇异点先看信号在哪一刻变了把压缩包名称拆开“qiyidian”是“奇异点”的拼音“小波模极大值”才是核心算法信号在某一瞬间发生跳变、脉冲或斜率突变时普通差分和FFT只能告诉你“有大变化”却回答不了变化发生在第几个采样点、突变有多尖锐。小波模极大值方法把原始信号放到多尺度空间里观察每个突变点的小波系数模极大值会沿尺度方向形成一条极值线其衰减速率由Lipschitz指数决定。把这个指数估算出来就能把“突变”从定性变成定量。下面这套流程面向故障诊断、信号处理和MATLAB工程实践从算法原理、代码实现、参数调整到验证技巧完整走一遍从原始波形到奇异点列表的路线。新手能按步骤复现熟手也能看到噪声条件下参数取舍的边界。2. 小波模极大值的数学基础Lipschitz指数与极值线传播2.1 Lipschitz指数度量奇异强度奇异点不是非黑即白的“跳变”同样的电压突降电网上相邻节点捕捉到的陡峭程度不同同样的边缘像素经图像退化后梯度也会变缓。Lipschitz指数把“陡峭程度”收敛成一个实数。若x(t)在t0附近满足|x(t) - x(t0)| ≤ K·|t - t0|^α则称x(t)在t0处具有Lipschitz指数α。α1对应光滑可微点α0对应有界间断阶跃α-1对应比阶跃更尖锐的冲击理想脉冲0和1之间的值对应斜坡状突变。指数越小奇异越强。在离散信号里直接算这个指数很困难因为采样点落在突变处的相位会严重影响差分结果。小波变换的优势在于它相当于用可伸缩的窗口对信号做带通卷积每个尺度的输出都携带该尺度下突变轮廓的信息。Mallat证明了小波系数模极大值满足log|Wx(s,t)| ≤ log K α·log s这给了工程上一个可操作的估计办法沿极值线取多个尺度的小波系数模在双对数坐标里拟合直线斜率即α。实现这一步所需的全部工作就是先把小波变换系数算对再把极值点从二维矩阵里准确拎出来。2.2 小波基的选择影响极值线含义先将数学结论落到小波基选择上。用于奇异点检测的小波必须满足两个条件具有一阶消失矩或以上能够反映信号的局部变化尽可能平滑避免小波自身的高频抖动在模极大值提取阶段产生伪点。最常用的是高斯函数的一阶导数即gauss1小波。gauss1在时域是单峰形态的疏波奇数对称模板本身没有直流分量卷积结果近似于对信号做一阶微分后再平滑。不同小波基对应不同的“极值位置含义”。表2-1列出的定位关系决定后续追踪逻辑建议在写代码前先想清楚自己的目标是找信号本身突变用gauss1还是找信号拐点用gauss2。表2-1 不同小波基对奇异点检测的定位差异小波基模极大值对应定位精度典型用途haar小波差分突变点一般快速粗检测gauss1信号一阶导数过零点高奇异点检测首选gauss2信号拐点中图像边缘检测morse时频脊线点低时频分析不推荐做奇异性估计为什么不建议直接用MATLAB的cwt默认morse小波估Lipschitz指数morse小波在尺度频带上的定位更偏向时频分析其相位特性和极值衰减规律比gauss1复杂拟合出的斜率波动大解释起来也更困难。gauss1的另一个好处是模板解析式简单便于在调试时手算特定尺度下的期望模值。如果需要确认当前MATLAB版本是否支持直接调用gauss1在命令行执行waveinfo(gaus)查看即可。2.3 极值线在尺度空间中的传播与偏移一个小窍门若在最小尺度上检测到模极大值并且随着尺度增大该极大值点持续存在、位置漂移不超过各自尺度下小波支撑半径则基本可判定这是一个真实奇异点。理论来源是奇异点周围小波系数构成锥形影响区极值线穿过该区域一直延伸到最大尺度。反过来噪声形成的极值线覆盖尺度很短通常在2到3个尺度就消失。定位和偏移的量级需要估计。尺度s下gauss1的等效支撑半径大约为2.5s个采样间隔若该尺度下有另一个同极性奇异点在附近两条极值线会在大尺度方向汇合。实际检测时时间坐标取极值线最细尺度端的位置而不是把整条线平均这样能避免大尺度偏移引入的定位误差。下面这段代码是极值线追踪和Lipschitz拟合的基础也解释了为什么会用双对数坐标scales 2 .^ (1:0.5:7); % 指数增长的尺度序列 amp 0.8 * scales .^ (-0.35); % 假设α-0.35 loglog(scales, amp, o-); xlabel(尺度 s); ylabel(模极大值 |W|);若模值随尺度呈幂律衰减在双对数坐标里就是一条直线直线斜率即α。写代码之前先造出这样一组数据做基准实现阶段会少走很多弯路。3. MATLAB实现从连续小波变换到模极大值提取3.1 构造测试信号与gauss1小波核先搭建一个最小可运行流程。测试信号应包含至少两种尺度的突变阶跃和冲击。阶跃突变对应α0冲击对应α≈-1二者在模极大值提取中的行为差异足够明显。Fs 1000; t (0:999) / Fs; x sin(2*pi*50*t); % 干净背景 x(300:end) x(300:end) 0.6; % 阶跃突变位于第300点 x(700) x(700) 1.0; % 冲击突变位于第700点信号构造本身很简单第300点开始整体抬高形成阶跃第700点是单点脉冲。测试信号的作用是给整个流程提供一个可对照的真值检测完成后能直接拿检测位置与构造位置做差定位误差一目了然。接下来生成gauss1小波核。常见做法是先固定时间轴分辨率dt再按尺度s决定核长度。这样生成的每个尺度核都覆盖了高斯函数的完整支撑范围不会因为核太短把高频细节切掉dt 1 / Fs; nmax ceil(8 * s / dt); % 尺度s下核的半长度 nu -nmax:nmax; psi -(nu * dt / s) .* exp(-(nu * dt / s).^2 / 2); W conv(x, psi, same) * (dt / sqrt(s));这段代码的核心是把连续小波定义离散化小波模板ψ(u)的自变量u取为(n·dt)/s尺度越大核越长对应观察的波形窗口越宽。末尾乘上dt/√s是把积分写成黎曼和后的小波能量归一化保证不同尺度下的模值处于同一量级后面Lipschitz拟合才有意义。完整循环就是逐尺度做卷积把每一层的输出作为一行堆成W矩阵。信号两端用conv的same模式会默认补零制造边界假极值。建议先对x做镜像延拓xext wextend(sym, 2, x, floor(nmax/2)2); Wext conv(xext, psi, same); W Wext(floor(nmax/2)3 : end - floor(nmax/2)-2);延拓后出现假极值的位置会被裁剪掉这是数字实现里最容易忽略的一层也是“同样代码别人结果好、自己结果乱”的最常见原因。完整代码中需要把延拓放在尺度循环外一次完成避免每个尺度重复延拓造成边界不一致。3.2 在每个尺度行上提取模极大值候选点小波系数矩阵建好后奇异点检测的下一步是把每个尺度上|W|的局部极大值抽出来。这里要区分两个概念单尺度内的局部极大值和跨尺度同一条极值线上的点。前者用一维滑动比较或islocalmax即可threshold 0.1 * max(abs(W(:))); cand cell(length(scales), 1); for k 1:length(scales) w abs(W(k, :)); lm islocalmax(w, MinProminence, 0.05 * max(w)); lm lm (w threshold); pos find(lm); cand{k} [pos(:), w(lm(:))]; % [位置模值] endislocalmax的MinProminence参数值得展开说。没有这个参数时只要某点两侧比邻点大就会判定为极大值含噪信号里会出现大量假峰。MinProminence设为当前尺度最大模值的5%能滤掉明显的浅峰。threshold阈值再兜底一次把整体模值很低的候选点全部去掉。如果你的MATLAB版本较老没有islocalmax可以用diff(w)完成等价的局部极大值判定核心逻辑不变。cand{k}里存放的是二维数组第一列是时间位置第二列是对应的模值。后续极值线追踪完全依赖这两个字段所以提取阶段的过滤宁可偏严也不要偏松否则追踪代码会因为这些太密的候选点而连线错乱。3.3 跨尺度极值线追踪追踪的目标是把同一奇异点在不同尺度产生的模极大值连成线。常见做法是从最小尺度出发向大尺度找数最近邻当前尺度的候选点pos在上一尺度所有候选点中找时间差最小的点若差值小于允许范围就判定为同一条线。允许范围一般取0.2倍当前尺度对应的小波支撑宽度或干脆固定为3到5个采样点dtmax 0.2 * scales(k) * dt 2 * dt; ids zeros(size(cand{k}, 1), 1); for k 2:length(scales) c1 cand{k-1}; c2 cand{k}; for i 1:size(c1, 1) [mind, j] min(abs(c2(:, 1) - c1(i, 1))); if mind dtmax id assign_line(c1(i,:), c2(j,:)); % 把两点并到一条线 end end end此处的assign_line是你要自己维护的线id分配函数在正式工程里可以用图论里的连通分量来做避免多条线在赋值时互相覆盖。更简单的实现是把lines存成元胞数组每条线保存posSeq和ampSeq匹配成功就在线尾追加节点如果当前候选点与任何已有线都不匹配则新开一条线。追踪方向也可以从大尺度向小尺度反向进行但对噪声信号从小尺度出发更容易保持定位精度。追踪完成后一条有效的极值线至少要覆盖60%以上尺寸。这个比例是过滤噪声伪线的基本标准具体取值在第4章结合参数展开。3.4 用双对数拟合估计Lipschitz指数每条线的posSeq对应的尺度系数已知ampSeq由cand第二列得到。接下来只需要取对数后做线性拟合function alpha est_lipschitz(scales, amps, kRange) p polyfit(log(scales(kRange)), log(amps(kRange)), 1); alpha p(1); end调用例子如果某条线在尺度索引5到12之间均有数据就令kRange5:12得到该点的α估计值。为什么要舍弃两端的尺度最小尺度的模值容易受噪声扰动最大尺度的模值又可能被邻近奇异点通过极值线合并而干扰中间段更接近“只受该奇异点控制”的理论区域拟合稳定性最高。实际操作时通常做两次过滤第一次按覆盖率过滤线第二次按拟合残差过滤。若拟合残差的均方根超过0.3说明该线位置附近可能有两个互相影响的奇异点或者小波变换数值精度有问题需要回到第3.1节的延拓部分检查。拟合得到α后按表3-1做物理含义映射。表3-1 估计α与奇异类型的对应参考估计α范围奇异类型典型对象-1.2 ~ -0.8冲击局部放电脉冲、撞击信号-0.2 ~ 0.2阶跃相位切换、机械断裂0.3 ~ 0.7斜坡切变趋势拐点、缓慢饱和0.8 以上近似光滑一般不作为奇异点处理4. 参数选择与噪声场景下的小波模极大值调节4.1 尺度范围与尺度步长怎么定尺度范围的选择同时影响计算量和检测精度。最小尺度决定定位精度尺度1对应的gauss1核只覆盖约8个采样点定位误差也在1-2个采样点之内。如果信号采样率是1000Hz最小尺度取1就够若信号本身带宽有限例如传感器前端内置低通滤波则需要把最小尺度抬高到2-4避免检测出大量高频伪突变。最大尺度随着噪声水平和信号长度的变化而变化。理论上尺度取到信号长度的1/4以上没有意义因为此时整个信号基本都在小波支撑内奇异点已无法在时域上区分。经验上取信号长度的1/16到1/8。表4-1给出典型的初始参数组合实际使用时按信号能量分布微调最大尺度即可。表4-1 不同采样率下的尺度范围建议Fs信号长度最小尺度最大尺度步长128Hz2s1160.51kHz1s1640.510kHz0.5s21280.550kHz0.1s4641步长决定极值线匹配成功率。指数增长序列2^(a:b:c)中b取0.5时相邻尺度之间小波核长度相差约1.4倍核形状变化平缓追踪稳定。b取1时计算量小但相邻尺度间的模值可能突变对噪声敏感。推荐先从b0.5开始若信号很短或者只需要做粗检测再换到b1。4.2 阈值与显著性过滤固定百分比threshold0.1*max(|W|)在背景噪声很小时工作正常一旦噪声水平升高真实微弱奇异点的模值可能远低于背景峰值固定百分比就直接漏检。更稳健的做法是从中位数估计噪声能量noiseLevel median(abs(W(:))) / 0.6745; thr 4 * noiseLevel;0.6745是高斯噪声标准差与中位数绝对偏差的换算常数4倍噪声标准差对应的虚警率已经很低。该方法不需要知道信号的真实信噪比对工程应用更友好。如果信号中存在明显的有色噪声中位数会被低频成分抬高推荐改用每一尺度行单独计算中位数得到每层不同的阈值。噪声本身也会在尺度平面里形成极值点但这些极值点的显著性随尺度增大迅速衰退。若某条极值线覆盖尺度数少于总数一半其拟合出的α可信度非常低建议直接丢弃。阈值和覆盖率两个参数要一起调节只动一个往往压不住噪声。4.3 伪极值线的判别覆盖率与线长追踪完成后按线覆盖的尺度数过滤valid false(length(lines), 1); for i 1:length(lines) valid(i) numel(lines(i).amp) 0.6 * length(scales); end lines lines(valid);0.6这个比例需要结合尺度范围做微调。当尺度步长b0.5、总尺度数是13时真实极值线通常能覆盖全部13层噪声线一般只覆盖2到4层。如果最小尺度或最大尺度设置过短真线的覆盖率也会下降因此判断标准以“至少覆盖总尺度数的60%”为准不要直接写死为某个固定数量。白噪声在尺度谱上还有一个明显特征极值线分布均匀且密度高。过滤时可以把尺度平面的二维密度也纳入考量利用histcounts2统计单位区域的线数密度超标的区域整块剔除。这样处理对突发性强噪声更有效因为突发噪声通常只集中在某几个尺度上。4.4 低信噪比下的预平滑替代方案当信噪比降到10dB以下时只靠提高阈值往往压制不住噪声极值线需要在对小波系数做极值提取之前先做时间方向的中值滤波。常见做法是对每一尺度行w做medfilt1(w, m)m取该尺度小波主瓣宽度的一半。中值滤波在保边缘的同时能显著降低噪声峰缺点是小尺度端的极值位置会偏移几个样本点定位精度会略微下降。另一种方案是提高阈值的同时放宽dtmax。噪声让极值线位置抖动如果追踪窗口过紧真实极值线会被切碎成若干短线。此时可以把dtmax从固定3个采样点调整为0.4*scales(k)dt3dt。放宽追踪窗口的代价是邻近两个奇异点的极值线可能被误并为一条解决方法是检查线内模值的平滑性若同一条线内模值出现非单调反复则拆线重连。5. 验证技巧用合成信号检验小波模极大值的整条链路5.1 构造三类标准奇异点写完代码后第一件事不是跑真实数据而是用理论α已知的三类信号验证。阶跃点的Lipschitz指数理论为0脉冲点理论接近-1斜坡突变的理论在0.5附近。构造这样一支信号后运行前述函数比对检测位置与估计α。t (0:999) / 1000; x zeros(1, 1000); x(400) 1; % 脉冲 x(600:end) x(600:end) 1; % 阶跃 x(800:end) x(800:end) (0:199) / 200; % 斜坡切变脉冲在单点注入能量阶跃从第600点开始持久改变电平斜坡则在第800点开始引入斜率。三者的理论α分别接近-1、0和0.5足够覆盖多数故障信号中的突变形态。构造完成后分别调用极值提取和追踪函数观察三条极值线是否都出现。5.2 与理论指标做定量对照验证时要看的三个指标定位误差、α估计误差、极值线覆盖率。定位误差应小于最小尺度下的核半径α误差在±0.15以内算通过覆盖率应接近100%。任何一个不合格都优先回溯到边界延拓与阈值设置而不是直接改拟合代码。把这三个指标打印出来作为回归基准之后每改动一次参数跑一遍合成信号对比结果能快速发现是参数还是实现出了问题。若定位误差持续偏大可以在极值线追踪后用一次一维抛物线插值对峰值位置做亚像素修正这是不增加计算量就能提升定位精度的常用手段。5.3 按qiyidian.rar的组织习惯拆分函数标题里的qiyidian.rar让人联想到一个现成的MATLAB包。不管解压后原始文件怎么命名按“数据读取-尺度变换-极值提取-线追踪-指数估计”拆成五个函数去做以后换数据、换小波基时只需改对应层。例如换用gauss2只需要替换尺度变换里的模板函数其他模块完全不动。推荐的调用顺序W wt_sgram(x, scales, Fs); cand wt_findmax(W, MinProminence, 0.05); lines wt_track(cand, scales, dtmax, 0.3); result wt_lipschitz(lines, scales, MinCoverage, 0.6);参数名在函数接口里显式暴露比在脚本里写散落的魔数更容易维护。把合成信号验证脚本放到最前面作为一个独立的testSmoke.m形成每次改代码后必跑一次的回归习惯。这样即使后面接入的工程信号形态复杂也能保证小波模极大值这层核心算法始终处于可信状态。本文还有配套的精品资源点击获取
