博客

目录

语法

说明

示例

寻找向量中的峰值

查找峰值及其位置

findpeaks函数的功能是寻找局部最大值。

语法

pks = findpeaks(data)
[pks,locs] = findpeaks(data)
[pks,locs,w,p] = findpeaks(data)
[___] = findpeaks(data,x)
[___] = findpeaks(data,Fs)
[___] = findpeaks(___,Name,Value)
findpeaks(___)

说明

pks = findpeaks(data)

返回带有输入信号向量data的局部最大值(峰值)的向量。局部峰值是大于其两个相邻样本或等于Inf的数据样本。峰值按出现的顺序输出。不包括非Inf信号端点。如果峰是平的，函数只返回指数最低的点。

[pks，locs] = findpeaks(data)

还会返回峰值出现的索引。

[pks，locs，w，p] = findpeaks(data)

还返回峰的宽度作为向量w，峰的突出度作为向量p。

[___] = findpeaks(data，x)

将x指定为位置向量，并返回以前语法中的任何输出参数。loc和w用x表示。

[___] = findpeaks(data，Fs)

指定数据的采样速率Fs。假设第一个数据样本是在时间零点采集的。loc和w转换为时间单位。

[___] = find peaks(_ _ _，Name，Value)

指定使用名称-值参数以及以前语法中的任何输入参数的选项。

不带输出参数的findpeaks(___)绘制信号并覆盖峰值。

示例

寻找向量中的峰值

定义一个有三个峰值的矢量并画出来。

data = [25 8 15 5 6 10 10 3 1 20 7];
plot(data)

如图所示：

求局部极大值。峰值按出现的顺序输出。第一个样本不包括在内，尽管它是最大值。对于平峰，该函数仅返回指数最低的点。

pks = findpeaks(data)


pks = 1×3

    15    10    20

使用不带输出参数的findpeaks来显示峰值。

findpeaks(data)

如图所示：

查找峰值及其位置

创建一个由一系列钟形曲线组成的信号。指定每条曲线的位置、高度和宽度。

x = linspace(0,1,1000);

Pos = [1 2 3 5 7 8]/10;
Hgt = [3 4 4 2 2 3];
Wdt = [2 6 3 3 4 6]/100;

for n = 1:length(Pos)
    Gauss(n,:) = Hgt(n)*exp(-((x - Pos(n))/Wdt(n)).^2);
end

PeakSig = sum(Gauss);

画出各条曲线及其总和。

plot(x,Gauss,'--',x,PeakSig)

如图所示：

使用默认设置的findpeaks查找信号的峰值及其位置。

[pks,locs] = findpeaks(PeakSig,x);

使用findpeaks绘制峰并标记它们。

findpeaks(PeakSig,x)

text(locs+.02,pks,num2str((1:numel(pks))'))

如图所示：

从最高到最低对山峰进行排序。

[psor,lsor] = findpeaks(PeakSig,x,'SortStr','descend');

findpeaks(PeakSig,x)

text(lsor+.02,psor,num2str((1:numel(psor))'))

如图所示：

峰值突起

创建一个信号，该信号由一个完整余弦周期上的钟形曲线之和组成。指定每条曲线的位置、高度和宽度。

x = linspace(0,1,1000);

base = 4*cos(2*pi*x);

Pos = [1 2 3 5 7 8]/10;
Hgt = [3 7 5 5 4 5];
Wdt = [1 3 3 4 2 3]/100;

for n = 1:length(Pos)
    Gauss(n,:) =  Hgt(n)*exp(-((x - Pos(n))/Wdt(n)).^2);
end

PeakSig = sum(Gauss)+base;

画出各条曲线及其总和。

plot(x,Gauss,'--',x,PeakSig,x,base)

如图所示：

使用findpeaks定位并绘制显著度至少为4的峰。

findpeaks(PeakSig,x,'MinPeakProminence',4,'Annotate','extents')

如图所示：

只有最高和最低的峰值满足该条件。 显示所有峰的日珥和半日珥宽度。

[pks,locs,widths,proms] = findpeaks(PeakSig,x);
widths
widths = 1×6

    0.0154    0.0431    0.0377    0.0625    0.0274    0.0409




proms
proms = 1×6

    2.6816    5.5773    3.1448    4.4171    2.9191    3.6363

寻找分离度最小的峰

太阳黑子是一种循环现象。已知它们的数量大约每11年达到一个峰值。 加载文件spocitizen . dat，其中包含了从1700年到1987年每年观测到的太阳黑子的平均数量。找到并画出最大值。

load sunspot.dat

year = sunspot(:,1);
avSpots = sunspot(:,2);

findpeaks(avSpots,year)

如图所示：

通过忽略彼此非常接近的峰值来改进对周期持续时间的估计。再次查找并绘制峰，但现在将可接受的峰-峰分离限制为大于6年的值。

findpeaks(avSpots,year,'MinPeakDistance',6)

如图所示：

使用findpeaks返回的峰位置来计算最大值之间的平均间隔。

[pks,locs] = findpeaks(avSpots,year,'MinPeakDistance',6);

meanCycle = mean(diff(locs))


meanCycle = 10.9600

使用年份数据创建日期时间数组。假设太阳黑子在每年的3月20日被计数，接近春分。找出太阳黑子的高峰期。使用年数函数将最小峰间距指定为持续时间。

ty = datetime(year,3,20);

[pk,lk] = findpeaks(avSpots,ty,'MinPeakDistance',years(6));

plot(ty,avSpots,lk,pk,'o')

如图所示：

使用日期时间函数计算平均太阳黑子周期。

dttmCycle = years(mean(diff(lk)))

dttmCycle = 10.9600

用这些数据创建一个时间表。以年为单位指定时间变量。绘制数据。显示时间表的最后五个条目。

TT = timetable(years(year),avSpots);
plot(TT.Time,TT.Variables)

如图所示：

entries = TT(end-4:end,:)


entries=5×1 timetable
      Time      avSpots
    ________    _______

    1983 yrs     66.6  
    1984 yrs     45.9  
    1985 yrs     17.9  
    1986 yrs     13.4  
    1987 yrs     29.3  

约束峰值特征

加载以7418 Hz采样的音频信号。选择200个样本。

load mtlb
select = mtlb(1001:1200);

找出间隔至少5 ms的峰。 要应用此约束，findpeaks会选择信号中最高的峰值，并消除其5毫秒内的所有峰值。然后，该函数对最高的剩余峰重复该过程，并进行迭代，直到用完要考虑的峰。

findpeaks(select,Fs,'MinPeakDistance',0.005)

如图所示：

找出幅度至少为1 V的峰值。

findpeaks(select,Fs,'MinPeakHeight',1)

如图所示：

找出比相邻样本至少高1 V的峰值。

findpeaks(select,Fs,'Threshold',1)

如图所示：

在信号达到更高值之前，找出任一侧下降至少1 V的峰值。

findpeaks(select,Fs,'MinPeakProminence',1)

如图所示：

饱和信号的峰值

如果数据大于给定的饱和点，传感器可以返回削波读数。可以选择忽略这些没有意义的峰值，或者将它们纳入您的分析中。 产生一个信号，该信号由嵌入方差为0.1的高斯白噪声中的频率为5 Hz和3 Hz的三角函数的乘积组成。信号以100 Hz的速率采样一秒钟。重置随机数生成器以获得可重现的结果。

rng default

fs = 1e2;
t = 0:1/fs:1-1/fs;

s = sin(2*pi*5*t).*sin(2*pi*3*t)+randn(size(t))/10;

通过截断每个大于指定界限0.32的读数来模拟饱和测量。绘制饱和信号。

bnd = 0.32;
s(s>bnd) = bnd;

plot(t,s)
xlabel('Time (s)')

如图所示：

找到信号的峰值。findpeaks只报告每个平峰的上升沿。

[pk,lc] = findpeaks(s,t);

hold on
plot(lc,pk,'x')

如图所示：

使用“阈值”名称-值对来排除平峰。要求最小幅度差为10^-4 一座山峰和它的邻居之间。

[pkt,lct] = findpeaks(s,t,'Threshold',1e-4);

plot(lct,pkt,'o','MarkerSize',12)

如图所示：

确定峰宽

创建一个由一系列钟形曲线组成的信号。指定每条曲线的位置、高度和宽度。

x = linspace(0,1,1000);

Pos = [1 2 3 5 7 8]/10;
Hgt = [4 4 2 2 2 3];
Wdt = [3 8 4 3 4 6]/100;

for n = 1:length(Pos)
    Gauss(n,:) =  Hgt(n)*exp(-((x - Pos(n))/Wdt(n)).^2);
end

PeakSig = sum(Gauss);

画出各条曲线及其总和。

plot(x,Gauss,'--',x,PeakSig)
grid

如图所示：

使用半日珥作为参考测量峰的宽度。

findpeaks(PeakSig,x,'Annotate','extents')

如图所示：

再次测量宽度，这次使用一半高度作为参考。

findpeaks(PeakSig,x,'Annotate','extents','WidthReference','halfheight')
title('Signal Peak Widths')

如图所示：