我试图找到一些振荡的阻尼率，我从一个模拟输出。
为此，我首先找到这些振荡的峰值，然后用指数拟合。
该方法适用于大多数情况，但对于我昨天遇到的一个情况，我发现该方法不能正常工作，并且没有明显的原因可以找到它不能工作。似乎有一些问题与我使用的峰值查找，不知何故错过了一些情况下的峰值。
下面是查找峰的程序：（注：列表Tarr、Zarr、Alpharr中的变量取自之前运行的细胞）

from scipy.optimize import curve_fit

def exp(t, A, lbda):
  return A * np.exp(-lbda * t)

def find_peaks(x, y):
    peak_x = []
    peak_vals = []
    for i in range(len(y)):
        if i == 0:
          if y[i+1]<y[i]:
            peak_x.append(x[i])
            peak_vals.append(y[i])
          else:
            continue
        if i == len(y) - 1:
            continue
        if (y[i-1] < y[i]) and (y[i+1] < y[i]):
            peak_x.append(x[i])
            peak_vals.append(y[i])
    return np.array(peak_x), np.array(peak_vals)

Tarr=[Tplot0,Tplot1,Tplot2,Tplot3,Tplot4]
Zarr=[Zrel0,Zrel1,Zrel2,Zrel3,Zrel4]
Alpharr=[alpha0,alpha1,alpha2,alpha3,alpha4]

for x0,y0,alf in zip(Tarr,Zarr,Alpharr):

  peak_x0,peak_y0=find_peaks(x0,y0)
  popt, pcov = curve_fit(exp, peak_x0, peak_y0)
  print(*[f"{val:.2f} +/- {err:.2f} __ " for val, err in zip(popt, np.sqrt(np.diag(pcov)))])
  print("decay rate = ",popt[1])
  plt.plot(x0,y0)
  plt.plot(peak_x0,exp(peak_x0,*popt),'-r', label=str((popt[0]))+"*exp(-"+str((popt[1]))+")")
  plt.xlabel("time")
  plt.ylabel("Zrel ")
  plt.title("n = "+str(n)+" ; alpha = "+str(alf))
  plt.legend()
  plt.show()

对比以下两种情况：

1）

2）

在第二种情况下，振荡在衰减到零之前增加，而在第一种情况下，这些振荡非常小（大约为2）。我本来希望我的算法几乎连接第二个图中的峰值（以给予最佳拟合）。但似乎由于某种原因，它错过了前两个峰值。