如何获得信号的高包络和低包络

问题描述

我有一个相当嘈杂的数据，我正在试图计算出信号的高包络和低包络。它有点像MATLAB中的这个例子：

http://uk.mathworks.com/help/signal/examples/signal-smoothing.html

在&Q；中提取峰值包络&Q；。Python中有没有类似的函数可以做到这一点呢？我的整个项目都是用Python编写的，最坏的情况是，我可以提取Numpy数组并将其放入MATLAB中，然后使用该示例。但我更喜欢马特普利布的样子。CBA在MATLAB和Python之间执行所有这些I/O.

谢谢，

推荐答案

据我所知，Numpy/Scipy/Python中没有这样的函数。然而，创建一个并不是那么困难。总体思路如下：

给定值的向量：

1. 查找(%s)的峰位置。我们称它们为(U)
2. 找到%s的槽的位置。我们称它们为(L)。
3. 将模型匹配到(U)值对。我们称它为(U_P)
4. 将模型与(L)值对相匹配。我们称它为(L_P)
5. 在(S)的区域上计算(U_P)以得到上包络的插值值。(我们称它们为(Q_U))
6. 在(S)的定义域上计算(L_P)以得到下包络的插值值。(我们称它们为(Q_L))。

如您所见，它由三个步骤组成(查找位置、拟合模型、评估模型)，但应用了两次，一次用于信封的上半部分，一次用于下半部分。

要收集(S)的"峰"，您需要定位(S)的斜率由正变负的点；要收集(S)的"谷"，您需要定位(S)的斜率由负变正的点。

峰值示例：S=[4，5，4]5-4为正4-5为负

一个低谷示例：S=[5，4，5]4-5为负5-4为正

以下是一个示例脚本，可帮助您开始使用大量的内联注释：

from numpy import array, sign, zeros
from scipy.interpolate import interp1d
from matplotlib.pyplot import plot,show,hold,grid

s = array([1,4,3,5,3,2,4,3,4,5,4,3,2,5,6,7,8,7,8]) #This is your noisy vector of values.

q_u = zeros(s.shape)
q_l = zeros(s.shape)

#Prepend the first value of (s) to the interpolating values. This forces the model to use the same starting point for both the upper and lower envelope models.

u_x = [0,]
u_y = [s[0],]

l_x = [0,]
l_y = [s[0],]

#Detect peaks and troughs and mark their location in u_x,u_y,l_x,l_y respectively.

for k in xrange(1,len(s)-1):
    if (sign(s[k]-s[k-1])==1) and (sign(s[k]-s[k+1])==1):
        u_x.append(k)
        u_y.append(s[k])

    if (sign(s[k]-s[k-1])==-1) and ((sign(s[k]-s[k+1]))==-1):
        l_x.append(k)
        l_y.append(s[k])

#Append the last value of (s) to the interpolating values. This forces the model to use the same ending point for both the upper and lower envelope models.

u_x.append(len(s)-1)
u_y.append(s[-1])

l_x.append(len(s)-1)
l_y.append(s[-1])

#Fit suitable models to the data. Here I am using cubic splines, similarly to the MATLAB example given in the question.

u_p = interp1d(u_x,u_y, kind = 'cubic',bounds_error = False, fill_value=0.0)
l_p = interp1d(l_x,l_y,kind = 'cubic',bounds_error = False, fill_value=0.0)

#Evaluate each model over the domain of (s)
for k in xrange(0,len(s)):
    q_u[k] = u_p(k)
    q_l[k] = l_p(k)

#Plot everything
plot(s);hold(True);plot(q_u,'r');plot(q_l,'g');grid(True);show()

这将产生以下输出：

需要进一步改进的地方：

1. 上述代码不会过滤可能出现在比某个阈值距离(TL)(例如时间)更近的波峰或波谷。这类似于envelope的第二个参数。不过，通过检查u_x,u_y的连续值之间的差异，添加它很容易。

2. 但是，比前面提到的更快的改进是，在内插上、下包络函数之前，使用移动平均过滤器对数据进行低通滤波。您可以很容易地通过卷积您的(S)与适当的移动平均过滤器。这里不需要详细介绍(如果需要，也可以这样做)，要生成一个在N个连续样本上运行的移动平均过滤器，您可以这样做：s_filtered = numpy.convolve(s, numpy.ones((1,N))/float(N)。(N)越高，您的数据显示得越平滑。但请注意，这将使您的(S)值(N/2)采样向右移动(在s_filtered中)，这是因为平滑过滤器的名称为group delay。有关移动平均线的更多信息，请参见this link。

希望这能有所帮助。

(如果提供了有关原始应用程序的更多信息，我很乐意更正您的答复。也许可以以更合适的方式对数据进行预处理(？)

这篇关于如何获得信号的高包络和低包络的文章就介绍到这了，希望我们推荐的答案对大家有所帮助，也希望大家多多支持！