请选择 进入手机版 | 继续访问电脑版

QQ登录

只需要一步,快速开始

 注册地址  找回密码
查看: 421|回复: 0

一种峰值检测算法——AMPD算法(C语言实现)

[复制链接]
字体大小: 正常 放大

584

主题

1

听众

1422

积分

该用户从未签到

发表于 2023-11-27 12:11 |显示全部楼层
|招呼Ta 关注Ta
一种峰值检测算法——AMPD算法(C语言实现)
本文算法的原始论文出处:Algorithms | Free Full-Text | An Efficient Algorithm for Automatic Peak Detection in Noisy Periodic and Quasi-Periodic Signals | HTML (mdpi.com)
有位老哥在知乎写了Python代码:python代码

在数字信号处理中,经常涉及到波峰查找算法,如振动信号分析,样条插值法求包络等。对于周期信号或者准周期信号,文章介绍了一种名为Automatic multiscale-based peak detection (AMPD),即自动多尺度峰值查找算法。
同时对非周期信号的效果也很nice,强烈安利!
其优势是:
(1)算法本身(几乎)对信号具有良好的自适应性,唯一的假设是信号是周期的或者准周期的;
(2)抗噪能力强,后面可以看到,对周期性的要求也不是很高。
原理不多讲,可以直接来看原文,也比较简单,就是用一个多尺度的滑动窗口去两侧进行比较,寻找局部最大值。

下面为该算法的C语言实现:
  1. // 寻找数组最小值的下标
  2. int argmin(int* index, int index_len)
  3. {
  4.   int min_index = 0;
  5.   int min = index[0];
  6.   for (int i = 1; i < index_len; i++)
  7.   {
  8.     if (index[i] < min)
  9.     {
  10.       min = index[i];
  11.       min_index = i;
  12.     }
  13.   }
  14.   return min_index;
  15. }

  16. //寻找极值点函数
  17. // data是存放数据的数组
  18. //index是存放峰值点下标的数组
  19. //len_index是峰值个数,即index数组长度
  20. void AMPD(double* data,int* index,int *len_index)
  21. {
  22.   int* p_data = (int*)malloc(sizeof(int) * size); //size可以最大为数组长度
  23.   int* arr_rowsum = (int*)malloc(sizeof(int) * size);
  24.   int min_index, max_window_length;
  25.   for (int i = 0; i < size; i++)
  26.   {
  27.     p_data[i] = 0;
  28.   }
  29.   for (int k = 1; k <= size / 2 + 1; k++)
  30.   {
  31.     int row_sum = 0;
  32.     for (int i = k; i <= size - k; i++)
  33.     {
  34.       if ((data[i] > data[i - k]) && (data[i] > data[i + k]))
  35.         row_sum -= 1;
  36.     }
  37.     *(arr_rowsum + k - 1) = row_sum;
  38.   }
  39.   /*for (int i = 0; i < size/2; i++)
  40.   {
  41.     printf("%d\n", arr_rowsum[i]);
  42.   }*/
  43.   min_index = argmin(arr_rowsum, size/2); //此处为最大的窗口
  44.   //printf("%d\n", min_index);
  45.   max_window_length = min_index;
  46.   for (int k = 1; k < max_window_length + 1;k++)
  47.   {
  48.     for (int i = 1; i < size - k; i++)
  49.     {
  50.       if ((data[i] > data[i - k]) && (data[i] > data[i + k]))
  51.         p_data[i] += 1;
  52.     }
  53.   }
  54.   for (int i_find = 0; i_find < size; i_find++)
  55.   {
  56.     if (p_data[i_find] == max_window_length)
  57.     {
  58.       index[*len_index] = i_find;
  59.       (*len_index) += 1;
  60.     }
  61.   }
  62.   free(p_data);
  63.   free(arr_rowsum);
  64. }
复制代码
周期信号效果原始文章内太多了,这里就不展示了。
下面来展示一波我自己使用的,非周期信号的峰值寻找情况

VeryCapture_20231127120013.jpg

可以根据AMPD后得到的峰值数组,进行三次样条插值进行原信号包络的获取。(参考MATLAB内envelope函数)

寻找波谷的话直接将原始数据翻转一下,就可以得到波谷的下标了。



zan
您需要登录后才可以回帖 登录 | 注册地址

qq
收缩
  • 电话咨询

  • 04714969085
fastpost

关于我们| 联系我们| 诚征英才| 对外合作| 产品服务| QQ

手机版|Archiver| |繁體中文 手机客户端  

蒙公网安备 15010502000194号

Powered by Discuz! X2.5   © 2001-2013 数学建模网-数学中国 ( 蒙ICP备14002410号-3 蒙BBS备-0002号 )     论坛法律顾问:王兆丰

GMT+8, 2024-2-22 10:05 , Processed in 0.430196 second(s), 54 queries .

回顶部