测量信号的峰值检测

我们使用数据采集卡从设备读取信号,将其信号增加到峰值,然后回落到接近原始值。 为了find峰值,我们目前在arrays中search最高读数,并使用索引来确定在我们的计算中使用的峰值的时间。

如果最高值是我们正在寻找的峰值,那么这个效果很好,但是如果设备工作不正常,我们可以看到第二个峰值可能高于初始峰值。 我们在90秒内从16台设备每秒读取10个读数。

我最初的想法是循环阅读检查,看看前一个和下一个点是否小于当前find一个峰值,并build立一个峰值arrays。 也许我们应该在当前位置的任何一边平均观察一些点,以便在系统中产生噪声。 这是继续进行还是有更好的技术?


我们使用LabVIEW,并且检查了LAVA论坛,并且有许多有趣的例子。 这是我们的testing软件的一部分,我们正在努力避免使用太多的非标准VI库,所以我希望对所涉及的过程/algorithm提供反馈,而不是特定的代码。

你可以尝试信号平均,即对于每个点,平均值与周围3个或更多点。 如果噪音很大,那么即使这样做也可能无济于事。

我意识到这是语言不可知的,但是猜测你正在使用LabView,LabView中有很多预先打包的信号处理VI,可以用来平滑和降低噪音。 NI论坛是获得更多专门帮助的好地方。

有很多很多经典的峰值检测方法,其中任何一个都可能工作。 你将不得不看到什么,特别是数据的质量。 这里是基本的描述:

  1. (x(0), y(0))(x(n), y(n))之间的任意两点之间加上y(i + 1) - y(i) 0 <= i < n并且将这个T (“旅行”)称为合适的小k ,并将R (“上升”)设置为y(n) - y(0) + k kT/R > 1表示峰值。 如果由于噪音引起的大行程不太可能,或者如果噪音围绕基曲线形状对称分布,则这可以正常工作。 对于您的应用程序,请以高于给定阈值的分数接受最早的峰值,或者分析每个上升值的行程曲线以获取更多有趣的属性。

  2. 使用匹配的filter来评分与标准峰形状的相似度(实质上,使用针对某种形状的归一化点积来获得相似度的余弦度量)

  3. 去卷积一个标准的峰值形状,并检查高值(虽然我经常发现2对噪音对于简单的仪器输出不那么敏感)。

  4. D((x0, y0), (x1, y1)) + D((x1, y1), (x2, y2)) > D((x0, y0),(x2, y2)) ,如果x0 < x1 < x2, y1 > 0.5 * (y0 + y2) ,检查欧几里得距离D((x0, y0), (x1, y1)) + D((x1, y1), (x2, y2)) > D((x0, y0),(x2, y2)) ,依赖于三angular不等式。 使用简单的比率将再次为您提供评分机制。

  5. 为你的数据拟合一个非常简单的2-高斯混合模型(例如,Numerical Recipes有一个很好的现成的代码块)。 走较早的高峰。 这将正确处理重叠的峰值。

  6. 在数据中find最简单的高斯曲线,柯西曲线,泊松曲线或者什么曲线。 在广泛的范围内评估此曲线,并在logging峰值位置之后从数据副本中减去此曲线。 重复。 以模型参数(标准偏差可能,但有些应用可能关心峰度或其他特征)的最早峰值符合一些标准。 当从数据中减去峰值时,留意遗留的文物。 最佳比赛可能由上面#2中提到的比赛得分决定。

我已经完成了你以前的工作:在DNA序列数据中find峰值,从测量曲线中估计出衍生物的峰值,并在直方图中find峰值。

我鼓励你仔细参加正确的基线。 维纳滤波或其他滤波或简单的直方图分析通常是在噪声存在的情况下进行基线的简单方法。

最后,如果你的数据通常是嘈杂的,而且你从数据卡中得到的数据是未引用的单端输出(甚至是引用的,而不是差分),并且如果你将大量的观察数据平均到每个数据点,观察并丢弃第一个和最后一个四分位,并平均剩下的。 有一些这样的exception消除策略可以是非常有用的。

这个问题已经被详细研究了。

TSpectrum*类ROOT (核/粒子物理分析工具)中有一套非常先进的实现。 该代码适用于一维到三维数据。

ROOT源代码是可用的,所以你可以抓住这个实现,如果你想。

从TSpectrum类文档:

本课程中使用的algorithm已经在以下参考文献中发表:

[1] M.Morhac等:背景消除方法的多维重合伽马射线谱。 Nuclear Instruments and Methods in Physics Research A 401(1997)113-132。

[2] M.Morhac等:有效的一维和二维金解卷积及其在伽马射线光谱分解中的应用。 Nuclear Instruments and Methods in Physics Research A 401(1997)385-408。

[3] M.Morhac等人:在多维重合伽玛射线谱中的峰的鉴定。 Nuclear Instruments and Methods in Research Physics A 443(2000),108-125。

这些文件是从没有NIM在线订阅的那些人的类文件链接的。


所做的简短版本是将直方图展平以消除噪声,然后通过平坦直方图中的蛮力检测局部最大值。

我想为这个线程贡献一个我自己开发的algorithm:

它基于分散原理:如果一个新的数据点与某个移动平均数有一个给定的x个标准偏差,则algorithm信号(也称为z分数 )。 该algorithm是非常强大的,因为它构造了一个单独的移动均值和偏差,使信号不会破坏门槛。 因此,将来的信号具有大致相同的精度,而不pipe之前的信号量如何。 该algorithm需要3个input: lag = the lag of the moving window ,门threshold = the z-score at which the algorithm signalsinfluence = the influence (between 0 and 1) of new signals on the mean and standard deviation threshold = the z-score at which the algorithm signals influence = the influence (between 0 and 1) of new signals on the mean and standard deviation 。 例如, lag 5将使用最后的5个观察来平滑数据。 如果数据点距移动平均值3.5标准偏差,则threshold 3.5将发出信号。 0.5的influence给信号一半的正常数据点的影响。 同样,0的influence完全忽略了重新计算新阈值的信号:因此0的影响是最可靠的选项。

它的工作原理如下:

伪代码

 # Let y be a vector of timeseries data of at least length lag+2 # Let mean() be a function that calculates the mean # Let std() be a function that calculates the standard deviaton # Let absolute() be the absolute value function # Settings (the ones below are examples: choose what is best for your data) set lag to 5; # lag 5 for the smoothing functions set threshold to 3.5; # 3.5 standard deviations for signal set influence to 0.5; # between 0 and 1, where 1 is normal influence, 0.5 is half # Initialise variables set signals to vector 0,...,0 of length of y; # Initialise signal results set filteredY to y(1,...,lag) # Initialise filtered series set avgFilter to null; # Initialise average filter set stdFilter to null; # Initialise std. filter set avgFilter(lag) to mean(y(1,...,lag)); # Initialise first value set stdFilter(lag) to std(y(1,...,lag)); # Initialise first value for i=lag+1,...,t do if absolute(y(i) - avgFilter(i-1)) > threshold*stdFilter(i-1) then if y(i) > avgFilter(i-1) set signals(i) to +1; # Positive signal else set signals(i) to -1; # Negative signal end # Adjust the filters set filteredY(i) to influence*y(i) + (1-influence)*filteredY(i-1); set avgFilter(i) to mean(filteredY(i-lag,i),lag); set stdFilter(i) to std(filteredY(i-lag,i),lag); else set signals(i) to 0; # No signal # Adjust the filters set filteredY(i) to y(i); set avgFilter(i) to mean(filteredY(i-lag,i),lag); set stdFilter(i) to std(filteredY(i-lag,i),lag); end end 

演示

演示鲁棒阈值算法

>原来的答案

这个方法基本上来自David Marr的书“Vision”

高斯模糊您的信号与您的峰值的预期的宽度。 这消除了噪音尖峰,你的阶段数据是完好无损的。

然后边缘检测(LOG会做)

那么你的边缘是特征的边缘(如峰)。 在峰的边缘之间查看,按大小对峰进行sorting,然后就完成了。

我已经使用这个变化,他们工作得很好。

我想你想把你的信号与一个预期的,典型的信号相互关联 。 但是,我研究信号处理已经有那么长时间了,但是我没有太在意。

你可以在你的逻辑中应用一些标准的Devision ,并注意到x%以上的峰值。

我不太了解仪器仪表,所以这可能是完全不切实际的,但是这又可能是一个有帮助的不同方向。 如果您知道读数如何失效,并且在出现这种故障的峰值之间存在一定的时间间隔,为什么不在每个时间间隔进行梯度下降。 如果下降将您带回您以前search过的地区,您可以放弃它。 根据采样表面的形状,这也可以帮助您find比search更快的峰值。

期望的峰值和不需要的第二峰值之间是否存在质的差异? 如果两个峰都是“尖锐的” – 即持续时间短 – 当在频域中查看信号(通过执行FFT),您将在大多数频带获得能量。 但是,如果“良好”峰值可靠地将能量存在于不存在于“不良”峰值中的频率,反之亦然,则可能能够以这种方式自动区分它们。