python - 如何使用Python区分信号的下降沿？

Question

我正在努力区分一些信号，以计算弹簧质量系统（地震仪）的自由周期振动和阻尼比。我正在使用Python作为主要处理程序。我需要做的是导入此信号，解析该信号并找到下降沿，然后从每个振荡的顶部返回一个峰列表作为列表，以便我可以计算阻尼比。一旦确定了数据集中振荡的位置，自由周期的计算就非常简单。

我最忙的是如何解析列表，识别下降沿，然后捕获每个Z0..Zn元素。一旦知道下降沿在哪里，就可以使用FFT相当容易地计算出振荡频率，但是如果我处理整个文件，在释放之前，系统通电会产生大量能量，有时可能会迫使算法启动代表接近DC偏移而不是实际振荡的超低频。 （在较高的阻尼比下，这可能只有四个或五个可测量的反弹，这尤其是一个问题）。

有没有人对我该如何做有一些想法？现在，下面的代码使用屏幕快照中信号的任意分配值。但是，我需要代码来计算这些值。另外，我还没有确定如何创建峰值列表以计算阻尼比h。非常欢迎您提供帮助以解决这些问题的想法。由于我在Stackoverflow中的声誉很高，因此在以下链接的示例屏幕快照中包含了我的信号：
（男孩，我希望这能奏效！）
第谷的采样信号->

https://github.com/tychoaussie/Sigcal_v1/blob/066faca7c3691af3f894310ffcf3bbb72d730601/Freeperiod_Damping%20Ratio.jpg

##########################################################
import os, sys, csv
from scipy import signal
from scipy.integrate import simps
import pylab as plt
import numpy as np
import scipy as sp
#
# code goes here that imports the data from the csv file into lists.
# One of those lists is called laser, a time-history list in counts from an analog-digital-converter
# sample rate is 130.28 samples / second.
#
                                      #            
                                      # Find the period of the observed signal
                                      #
delta = 0.00767                       # Calculated elsewhere in code, represents 130.28 samples/sec
                                      # laser is a list with about 20,000 elements representing time-history data from a laser position sensor
                                      # The release of the mass and system response occurs starting at sample number 2400 in this particular instance.

sense = signal.detrend(laser[2400:(2400+8192)]) # Remove the mean of the signal
N = len(sense)
W    = np.fft.fft(sense)
freq = np.fft.fftfreq(len(sense),delta) # First value represents the number of samples and delta is the sample rate

                                      #
                                      # Take the sample with the largest amplitude as our center frequency. 
                                      # This only works if the signal is heavily sinusoidal and stationary 
                                      # in nature, like our calibration data.
                                      #

idx = np.where(abs(W)==max(np.abs(W)))[0][-1]
Frequency = abs(freq[idx]) # Frequency in Hz
period = 1/(Frequency*delta) # represents the number of samples for one cycle of the test signal.
                                      #
                                      # create an axis representing time.
                                      #

dt = [] # Create an x axis that represents elapsed time in seconds. delta = seconds per sample, i represents sample count
for i in range(0,len(sensor)):
    dt.append(i*delta)


                                      #
                                      # At this point, we know the frequency interval, the delta, and we have the arrays 
                                      # for signal and laser. We can now discriminate out the peaks of each rebound and use them to process the damping ratio of either the 'undamped' system or the damping ratio of the 'electrically damped' system.
                                      #

print 'Frequency calcuated to ',Frequency,' Hz.'

Answer 1

这是一个有点不合常规的想法，我认为这可能相当健壮，并且不需要大量的试探法和猜测。您所拥有的数据确实是高质量的，并且符合已知曲线，因此在这里很有帮助。在这里，我假设曲线的“好部分”具有以下形式：

V = a * exp(-γ * t) * cos(2 * π * f * t + φ) + V0     # [Eq1]

V：电压
t：时间
γ：阻尼常数
f：频率
a：起始振幅
φ：开始阶段
V0：直流偏移


算法概述

摆脱偏移

首先， calculate the derivative numerically。由于数据质量很高，因此噪声不会对事物产生太大影响。

电压导数 V_deriv具有与原始数据相同的形式：相同的频率和阻尼常数，但是相位 ψ和幅度 b不同，

V_deriv = b * exp(-γ * t) * cos(2 * π * f * t + ψ)         # [Eq2]

令人高兴的是，这将自动消除直流偏移。

备选方案：不需要完全执行此步骤-偏移量是曲线拟合的一个相对较小的复杂问题，因为您始终可以通过对整个曲线进行平均来很好地估计出偏移量。权衡是，如果不使用导数，则将获得更好的信噪比。

初步猜测

现在考虑导数曲线（如果跳过了最后一步，则考虑原始曲线）。从最右边的数据点开始，然后沿曲线向左尽可能多地振荡，直到达到振幅大于某个振幅阈值 A的振荡为止。您想找到曲线的一部分，其中包含一个振荡，具有良好的信噪比。

如何确定振幅阈值很难说。这取决于您的传感器的性能。我建议将此作为参数进行调整。

曲线拟合

现在，您已经捕获了一个振荡，您可以做很多简单的事情：估算频率 f和阻尼常数 γ。使用这些初始估计，您可以在此振荡右侧（包括振荡本身）的所有数据执行 nonlinear curve fit。

如果使用导数，则适合的函数是[Eq2]，如果使用原始曲线，则适合的函数是[Eq1]（但是无论如何它们都是相同的函数，所以要讨论点）。

为什么需要这些初始估算？对于适合于衰减波的非线性曲线，至关重要的一点是您必须对参数（尤其是频率和阻尼常数）做出良好的初步猜测。其他参数相对来说不太重要（至少这是我的经验），但是这是为了防万一，您可以通过以下方法获取其他参数：


如果从最大值开始，则阶段为 0；如果从最小值开始，则阶段为 π。
您也可以猜测起始振幅，尽管即使您将其设置为 1，拟合算法通常也会很好。


寻找边缘

此时的曲线拟合应该非常好–您可以知道，因为 discrepancy between your curve and the fit非常低。现在的诀窍是尝试向左增加曲线的范围。

如果您停留在正弦曲线区域内，则曲线拟合应保持良好（如何判断“良好”需要做一些实验）。但是，一旦触及曲线的平坦区域，误差将开始急剧增加，并且参数将开始偏离。这可以用来确定“好”数据在哪里结束。

但是，您不必逐点执行此操作–效率很低。 binary search在这里应该可以很好地工作（可能是它的“单面”变化）。

摘要

您不必遵循此确切过程，但基本要点是，您可以从最右边开始对一小部分数据执行一些分析，然后逐渐增加时域，直到发现一个点为止。错误正在变得越来越大，而不是应该变得更小。

此外，您还可以组合不同的启发式方法，以查看它们是否彼此一致。如果没有，则可能是一个相当粗略的数据点，需要一些手动干预。

请注意，我上面概述的算法的一个特殊优点是，您将在过程中获得所需的结果（阻尼常数和频率）以及不确定性估计。

我忽略了提及大多数复杂的数学和算法细节，以便提供总体概述，但是如果需要，我可以提供更多细节。