第十二篇：智能峰值/谷值检测算法详解

在前几篇中，我们专注于语音信号的预处理和预测延拓，但实际应用中我们常需要从波形中自动定位极值点——例如基音周期中寻找波峰、频谱中识别谐波峰值、或心跳信号中检测 R 波峰值。MATLAB 自带的 findpeaks 函数（需 Signal Processing Toolbox）功能强大，但有时我们需要更轻量、更可控的实现。findpeakm.m 提供了完全自包含的峰值/谷值检测算法，支持二次插值亚样本精度、宽度容差去重，并且无需任何工具箱依赖。本篇将深入剖析其核心算法，从一阶差分定位到插值优化，再到局部峰值去重，带你全面掌握这个实用工具。

1. 峰值检测的核心思想

峰值检测的本质是在离散序列中找到局部极大值点。一个点 x(k)

是峰值，当且仅当其左侧相邻点小于或等于它，且右侧相邻点也小于或等于它（对于严格的峰值，两侧严格小于）。但实际信号常伴有噪声，直接比较容易出现虚假极值。

基本流程：

1. 计算相邻点的差分（一阶导数近似）；
2. 找出差分由正变负的位置 —— 这些位置就是峰值候选点；
3. 若信号包含谷值，则可将信号取反后再执行上述流程；
4. 可选：对候选点进行二次插值，获得亚样本级精度的位置和幅值；
5. 可选：根据距离容差剔除相距太近的峰值（保留幅值较大者）。

findpeakm.m 完整实现了上述逻辑，并支持 ‘q’ 插值模式和 ‘v’ 谷值模式。

2. 代码结构总览

function[k,v]=findpeakm(x,m,w)% 输入：x - 信号向量% m - 模式字符串：'q' 表示二次插值，'v' 表示寻找谷值% w - 宽度容差（样本数），若两个峰值距离 ≤ w，则剔除较低的% 输出：k - 峰值位置（若 'q' 模式则为浮点数）% v - 峰值幅值

主要步骤：

· 将信号转为列向量，若为谷值模式则取反；
· 计算差分 dx = x(2:end) - x(1:end-1)；
· 找到上升段和下降段的索引（dx>0 为上升，dx<0 为下降）；
· 利用上升/下降段的起始和结束位置，确定峰值所在位置（处理平坦区）；
· 若为 ‘q’ 模式，用二次插值修正位置和幅值；
· 若指定了 w，则进行邻近峰值去重；
· 若为谷值模式，将幅值取反回来；
· 若无输出参数，则自动绘图。

3. 差分法定位峰值——核心算法详解

3.1 计算差分并标记上升/下降

dx=x(2:end)-x(1:end-1);r=find(dx>0);% 上升段的起始索引（指向上方点的索引）f=find(dx<0);% 下降段的起始索引

dx(i) = x(i+1) - x(i)。若 dx(i) > 0，说明从 i 到 i+1 是上升趋势；若 dx(i) < 0，则为下降。

3.2 计算相对于上升和下降的“时间距离”

这段代码是 findpeakm 的精髓，它通过累计“自上次上升/下降以来的样本数”来定位峰值。

dr=r;dr(2:end)=r(2:end)-r(1:end-1);rc=repmat(1,nx,1);rc(r+1)=1-dr;rc(1)=0;rs=cumsum(rc);% rs 向量：每个样本点距离最近一次上升点的样本数

类似地，fs 计算距离最近一次下降点的样本数。

3.3 确定峰值候选位置

峰值应满足：

· 它离最近一次上升点很近（rs < fs）：说明这个点正处于上升之后；
· 它离最近一次下降点也很近（fq < rq）：说明它即将转为下降；
· 并且 floor((fq - rs)/2) == 0：这个条件确保了峰值位于平坦区的中心（若存在平坦段）。

k=find((rs<fs)&(fq<rq)&(floor((fq-rs)/2)==0));v=x(k);

在没有平坦区时，fq - rs 通常为 1，此时 floor(0.5)=0，满足条件。如果出现一个平台（plateau），例如 [1, 2, 2, 1]，该条件会将峰值定位在平台的中心位置（第 2 或第 3 个点，取决于 floor 的行为）。

4. 二次插值 —— 实现亚样本精度

当信号峰值不是恰好落在采样点上时，我们可以用抛物线拟合三个相邻点（左、候选、右）来精确估计峰值的真实位置和幅值。

设候选点索引为 k

，其幅值为 x_k

，左右点为 x_{k-1}

和 x_{k+1}

。构造二次多项式：

f(t) = a t^2 + b t + c

令 t=0

对应候选点 k

，则：

· f(0) = c = x_k
· f(-1) = a - b + c = x_{k-1}
· f(1) = a + b + c = x_{k+1}

解得：
a = \frac{x_{k-1} + x_{k+1}}{2} - x_k

b = \frac{x_{k+1} - x_{k-1}}{2}

极值点位置（相对于候选点）为 t_{\max} = -\frac{b}{2a}

，对应的幅值为：
f(t_{\max}) = x_k - \frac{b^2}{4a}

当 a > 0

时，抛物线开口向上，那是极小值，但峰值处应为 a < 0

（开口向下）。若 a \approx 0

，说明为平坦区，则取中心。

代码实现：

ifany(m=='q')b=0.5*(x(k+1)-x(k-1));a=x(k)-b-x(k-1);j=(a>0);% 通常 a<0，此处 j 用于区分平坦区v(j)=x(k(j))+0.25*b(j).^2./a(j);k(j)=k(j)+0.5*b(j)./a(j);k(~j)=k(~j)+(fq(k(~j))-rs(k(~j)))/2;% 平坦区取中心end

注意：这里 a > 0 是异常情况（极小值），但实际峰值处 a 应为负。代码用 a 判断是否为平坦区（a 接近 0），如果 a>0 则强制按极小值修正，但通常不会发生。更稳健的实现应检查 a < -eps 才进行插值。

5. 邻近峰值去重（宽度容差 w）

当两个峰值之间的距离小于等于 w 时，我们只保留幅值较大的那个。这是一个非极大值抑制（NMS）过程。

ifnargin>2j=find(k(2:end)-k(1:end-1)<=w);whileany(j)j=j+(v(j)>=v(j+1));% 若前一个更高，则删除后一个；否则删除前一个k(j)=[];v(j)=[];j=find(k(2:end)-k(1:end-1)<=w);endend

这段代码非常巧妙：j 是相邻峰值距离不足 w 的位置索引，然后根据幅值比较，决定删除前一个还是后一个（j 指向较低者）。循环直至所有距离均大于 w。

6. 谷值检测模式

若 m 中包含 ‘v’，则函数先将信号取反（x = -x），执行完上述峰值检测后，再将幅值取反回来（v = -v）。这样谷值就变成了峰值，复用了同一套逻辑。

7. 实战示例：检测正弦波中的峰值

t=0:0.01:1;x=sin(2*pi*5*t)+0.1*randn(size(t));[k,v]=findpeakm(x,'q',0.1);% 二次插值，去重宽度 0.1 样本% 绘制结果findpeakm(x,'q',0.1);% 无输出时自动绘图

你会看到峰值位置被精确标记，且幅值接近 1。二次插值使得位置误差远小于采样间隔。

8. 与其他函数的关系

· 在语音分析中，峰值检测可用于基音周期提取（检测相邻波峰距离）。
· 在频谱分析中，检测谐波峰值可进行共振峰估计。
· 与 findSegment 结合，可在每个有效语音段内单独检测峰值，避免静音段的伪峰。

9. 本讲小结

· 我们深入剖析了 findpeakm.m 的核心算法：差分定位、二次插值、邻近去重。
· 理解了谷值模式、平坦区处理、亚样本精度等高级特性。
· 通过示例展示了其简单易用的接口和强大的绘图功能。

峰值检测是信号分析中极为常见的基础操作，掌握这个工具将使你在处理各种波形特征提取时游刃有余。下一篇，我们将从峰值检测回到语音端点分割，学习如何将帧级 VAD 标记转换为连续的语音段结构体。

📥 所有代码均已打包，点击下方链接免费获取：
下载链接

下篇预告：语音端点检测后处理——连续有话段分割。我们将进入 findSegment.m，了解如何将离散的 0/1 标签聚合成有意义的语音片段，并计算每个片段的时长，为后续的语音识别或特征提取准备结构化数据。敬请期待！

思考题：如果信号中有一个很宽的平坦峰值（如方波的顶部），findpeakm 会定位在平台中心。若你想要检测平台两侧的边缘，应如何修改算法？欢迎评论区讨论。