简单 scipy curve_fit 测试未返回预期结果

Question

我试图根据仅几个周期的测量来估计大约 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 度的相移 window。

如有任何想法，我们将不胜感激！

Answer 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

简单 scipy curve_fit 测试未返回预期结果

Simple scipy curve_fit test not returning expected results

curve-fitting

scipy

frequency-analysis