【问题标题】:Simple scipy curve_fit test not returning expected results简单的 scipy curve_fit 测试未返回预期结果
【发布时间】:2019-02-26 03:44:05
【问题描述】:

我试图根据仅几个周期的测量来估计大约 50Hz 的输入信号的幅度、频率和相位。频率需要精确到 0.01 Hz。由于信号本身将是一个非常清晰的正弦波,我正在尝试使用 SciPy 的 curve_fit 进行参数拟合。之前没用过,所以写了一个快速测试函数。

我首先生成虚拟余弦波的单个周期的样本

from math import *
import numpy as np
from scipy.optimize import curve_fit
import matplotlib.pyplot as plt

fs = 1000 # Sampling rate (Hz)
T = .1 # Length of collection (s)
windowlength = int(fs*T) # Number of samples
f0 = 10 # Fundamental frequency of our wave (Hz)

wave = [0]*windowlength
for x in range(windowlength):
    wave[x] = cos(2*pi*f0*x/fs)

t = np.linspace(0,T,int(fs*T)) # This will be our x-axis for plotting

然后我尝试将这些示例拟合到一个函数中,改编 scipy 提供的官方示例中的代码:https://docs.scipy.org/doc/scipy/reference/generated/scipy.optimize.curve_fit.html

# Define function to fit
def sinefit(x, A, ph, f):
    return A * np.sin(2*pi*f * x + ph)

# Call curve_fit
popt,cov = curve_fit(sinefit, t, wave, p0=[1,np.pi/2,10])

# Plot the result
plt.plot(t, wave, 'b-', label='data')
plt.plot(t, sinefit(t, *popt), 'r-', label='fit')

print("[Amplitude,phase,frequency]")
print(popt)

这给了我 popt = [1., 1.57079633, 9.9] 和情节

plot output

我的问题是:为什么我的频率关闭了?我用余弦波的确切参数初始化了curve_fit函数,那么LM算法的第一次迭代不应该意识到残差为零并且它已经得到了正确的答案吗?幅度和相位似乎是这样,但频率 0.1Hz 太低了。

我认为这是一个愚蠢的编码错误,因为原始波浪和拟合在情节中明确排列。我还确认它们之间的差异在整个样本中为零。如果它们真的是 0.1 Hz 异相,那么在我的 100 毫秒窗口内会有 3.6 度的相移。

任何想法将不胜感激!

【问题讨论】:

    标签: scipy curve-fitting frequency-analysis


    【解决方案1】:

    问题是您的数组t 不正确。 t 中的最后一个值是 0.1,但采样周期为 1/fs = 0.001,t 中的最后一个值应该是 0.099。即100个样本的次数为[0, 0.001, 0.002, ..., 0.098, 0.099]。

    您可以使用任一正确创建t

    t = np.linspace(0, T, int(fs*T), endpoint=False)
    

    或

    t = np.arange(windowlength)/fs  # Use float(fs) if you are using Python 2
    

    【讨论】:

      猜你喜欢
      • 1970-01-01
      • 1970-01-01
      • 2016-08-05
      • 1970-01-01
      • 1970-01-01
      • 1970-01-01
      • 2011-08-09
      • 2010-09-22
      • 2021-05-17
      相关资源
      最近更新 更多