python scipy.signal.peak_widths --> 绝对高度?(fft -3dB 阻尼)

inc*_*ito 2 python scipy

https://docs.scipy.org/doc/scipy/reference/ generated/scipy.signal.peak_widths.html

我认为链接函数只能计算相对高度处的峰宽度。有谁知道是否有一个函数可以计算所有峰值的固定值(peak_amplitude - x)的宽度?

目前我正在尝试更改原始内部函数“_peak_widths”。cimport 已经失败。这里只理解部分源代码。我在代码中添加了我要进行修改的地方。

 with nogil:
    for p in range(peaks.shape[0]):
        i_min = left_bases[p]
        i_max = right_bases[p]
        peak = peaks[p]
        # Validate bounds and order
        if not 0 <= i_min <= peak <= i_max < x.shape[0]:
            with gil:
                raise ValueError("prominence data is invalid for peak {}"
                                 .format(peak))
        height = width_heights[p] = x[peak] - prominences[p] * rel_height 
Run Code Online (Sandbox Code Playgroud)

将此处更改为 x[peak] - 3

        # Find intersection point on left side
        i = peak
        while i_min < i and height < x[i]:
            i -= 1
        left_ip = <np.float64_t>i
        if x[i] < height:
            # Interpolate if true intersection height is between samples
            left_ip += (height - x[i]) / (x[i + 1] - x[i])

        # Find intersection point on right side
        i = peak
        while i < i_max and height < x[i]:
            i += 1
        right_ip = <np.float64_t>i
        if  x[i] < height:
            # Interpolate if true intersection height is between samples
            right_ip -= (height - x[i]) / (x[i - 1] - x[i])

        widths[p] = right_ip - left_ip
        if widths[p] == 0:
            show_warning = True
        left_ips[p] = left_ip
        right_ips[p] = right_ip
Run Code Online (Sandbox Code Playgroud)

lag*_*gru 5

如果这仍然与您相关,您可以“按原样”使用scipy.signal.peak_widths通过传入修改来实现您想要的prominence_data。根据您自己的回答

import numpy as np
from scipy.signal import find_peaks, peak_prominences, peak_widths

# Create sample data
x = np.linspace(0, 6 * np.pi, 1000)
x = np.sin(x) + 0.6 * np.sin(2.6 * x)

# Find peaks
peaks, _ = find_peaks(x)
prominences, left_bases, right_bases = peak_prominences(x, peaks)
Run Code Online (Sandbox Code Playgroud)

peak_widths如文档 中所述,测量宽度的高度计算如下h_eval = h_peak - prominence * relative_height

我们可以通过参数prominence_data和来控制后两个变量rel_height。因此,prominence我们可以创建一个所有值都相同的数组,并使用它来创建绝对高度,而不是传递每个峰值不同的计算值:

# Create constant offset as a replacement for prominences
offset = np.ones_like(prominences)

# Calculate widths at x[peaks] - offset * rel_height
widths, h_eval, left_ips, right_ips = peak_widths(
    x, peaks, 
    rel_height=1,
    prominence_data=(offset, left_bases, right_bases)
)

# Check that h_eval is 1 everywhere
np.testing.assert_equal(x[peaks] - h_eval, 1)

# Visualize result
import matplotlib.pyplot as plt
plt.plot(x)
plt.plot(peaks, x[peaks], "x")
plt.hlines(h_eval, left_ips, right_ips, color="C2")
plt.show()
Run Code Online (Sandbox Code Playgroud)

在此输入图像描述

正如您所看到的,每个峰的宽度均以相同的恒定偏移量 1 进行评估。通过使用原始数据left_bases和我们right_bases提供的peak_prominences限制最大测量宽度(例如,参见 299 和 533 处的峰)。如果您想消除该限制,您必须自己创建这些数组。