python scipy.signal.peak_widths -> абсолютная высота?(демпфирование -3 дБ)

Question

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 

ИЗМЕНИТЬ ЗДЕСЬ x [пик] - 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

Answer 1

Я только что удалил содержимое c.Вот мое решение:

def gauss(x, p): # p[0]==mean, p[1]==stdev
    return 1.0/(p[1]*np.sqrt(2*np.pi))*np.exp(-(x-p[0])**2/(2*p[1]**2))

def _peak_widths(x,peaks,prop,val=3):

    i_min = prop['left_bases']
    i_max = prop['right_bases']
    peak = peaks[0]
    # Validate bounds and order
    height = x[peak] - val

    # Find intersection point on left side
    i = peak
    while i_min < i and height < x[i]:
        i -= 1
    left_ip = 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 = i
    if  x[i] < height:
        # Interpolate if true intersection height is between samples
        right_ip -= (height - x[i]) / (x[i - 1] - x[i])

    widths = right_ip - left_ip
    left_ips = left_ip
    right_ips = right_ip

    return [height, widths, int(left_ips), int(right_ips)]

if __name__ == '__main__':

    # Create some sample data
    known_param = np.array([2.0, 0.07])
    xmin,xmax = -1.0, 5.0
    N = 1000
    X = np.linspace(xmin,xmax,N)
    Y = gauss(X, known_param)
    fig, ax= plt.subplots()
    ax.plot(X,Y)

    #find peaks
    peaks, prop = signal.find_peaks(Y, prominence = 3.1)
    ax.scatter(X[peaks],Y[peaks], color='r')

    #calculate peak width
    y, widths, x1, x2 = _peak_widths(Y,peaks, prop)

    print(f'width = { X[x1] - X[x2]}')

    l = mlines.Line2D([X[x1],X[x2]], [y,y], color='r')
    ax.add_line(l)
    plt.show()