Scipy curve_fit可变参数拟合m阶二维傅里叶级数报错解决

数据老陈

2026-8-26

问题背景

我正在开展基于二维傅里叶级数的优化问题研究，需要实现如下m阶二维傅里叶级数公式：

目标是实现阶数m可自定义的拟合函数：输入指定阶数m即可对数据完成m阶二维傅里叶级数拟合，求解得到除已知的C_00外的最优级数系数。

使用Scipy库的curve_fit函数拟合时存在核心障碍：模型参数数量随阶数m动态变化，例如m=1时需拟合8个级数系数+2个角速度参数，m=2时需拟合24个级数系数+2个角速度参数，固定参数列表的常规写法无法适配动态参数数量，先后编写三个版本的拟合函数均触发报错。

报错复现

第一版函数报错

第一版核心代码如下：

c00=50
def fourier(x, v, w1, w2):
   f=c00
   k=0
   for i in range(1, m+1):
       f=f+v[k]*np.cos(i*w1*x)+v[k+1]*np.sin(i*w1*x)
       k+=2
   for j in range(1, m+1):
      for i in range(-m, m+1):
           f=f+v[k]*np.cos(i*w1*x+j*w2*x)+v[k+1]*np.sin(i*w1*x+j*w2*x)
           k+=2
return f

time = [] #从txt文件读取的x轴数据
energy = [] # 从txt文件读取的y轴数据
popt, pcov = curve_fit(fourier, time, energy)

运行后触发报错：

File "/home/arcalinux/Desktop/fourier_series.py", line 63, in fourier
f=f+v[k]np.cos(iw1x)+v[k+1]np.sin(iw1x)
IndexError: invalid index to scalar variable.

第二版函数报错

第二版函数代码如下：

def fourier2(args):
   f=c00
   k=3
   for i in range(1, m+1):
       f=f+args[k]*np.cos(i*args[1]*args[0])+args[k+1]*np.sin(i*args[1]*args[0])
       k+=2
   for j in range(1, m+1):
       for i in range(-m, m+1):
          f=f+args[k]*np.cos(i*args[0]*x+j*args[1]*x)+\
                          +args[k+1]*np.sin(i*args[0]*x+j*args[1]*x)
          k*=2
   return f

调用curve_fit时触发报错：

ValueError: Unable to determine number of fit parameters.

手动传参验证核心计算逻辑可正常运行，例如传入测试参数：

v=[1,1,1,1,1,1,1,1,1,1,1,1,1,1,1,1,1,1,1,1,1,1,1,1]
m=2
x=2
w1=3
w2=4
print(fourier(x,v,w1,w2))

可正常输出计算结果：

49.136961626348295

fourier2函数手动传参也可正常输出结果，项目完整代码如下：

#导入依赖库
import numpy as np
import math
from scipy.optimize import curve_fit
from statistics import mean
import matplotlib.pyplot as plt


m=2 #傅里叶级数阶数
number_of_coeff=(2*m+1)**2 - 1 
N=number_of_coeff
#构建级数
coefficients=[]

#第一部分，j=0
for i in range(1,m+1):
   coefficients.append('c'+str(i)+'0')
   coefficients.append('s'+str(i)+'0')

#第二部分
for j in range(1, m+1):
   for i in range(-m, 0):
       coefficients.append('g'+str(-i)+str(j))
       coefficients.append('t'+str(-i)+str(j))
   for i in range(0, m+1):
       coefficients.append('c'+str(i)+str(j))
       coefficients.append('s'+str(i)+str(j))
print(len(coefficients))
v=np.ones(N)
m=2
x=2
w1=3
w2=4
arg=[2,3,4,1,1,1,1,1,1,1,1,1,1,1,1,1,1,1,1,1,1,1,1,1,1,1,1]
arg3=[3,4,1,1,1,1,1,1,1,1,1,1,1,1,1,1,1,1,1,1,1,1,1,1,1,1]

c00=50
def fourier(x, v, w1, w2):
    print(type(v))
    print('v is ',v)
    f=c00
    k=0
    print(v[0])
    for i in range(1, m+1):
        f=f+v[k]*np.cos(i*w1*x)+v[k+1]*np.sin(i*w1*x)
    
        k+=2
    for j in range(1, m+1):
        for i in range(-m, m+1):
            f=f+v[k]*np.cos(i*w1*x+j*w2*x)+\
                v[k+1]*np.sin(i*w1*x+j*w2*x)
            k+=2
return f
def fourier2(args):
    f=c00
    k=3
    for i in range(1, m+1):
        f=f+args[k]*np.cos(i*args[1]*args[0])+\
            args[k+1]*np.sin(i*args[1]*args[0])
        k+=2
    for j in range(1, m+1):
        for i in range(-m, m+1):
            f=f+args[k]*np.cos(i*args[1]*args[0]+\
                j*args[2]*args[0])+\
                args[k+1]*np.sin(i*args[1]*args[0]+\
                j*args[2]*args[0])
        k+=2
return f
def fourier3(x,args):
    f=c00
    k=2
    for i in range(1, m+1):
        f=f+args[k]*np.cos(i*args[0]*x)+\
            args[k+1]*np.sin(i*args[0]*x)
        k+=2
    for j in range(1, m+1):
        for i in range(-m, m+1):
            f=f+args[k]*np.cos(i*args[0]*x+j*args[1]*x)+\
                args[k+1]*np.sin(i*args[0]*x+j*args[1]*x)
            k+=2
return f

#读取数据
time = [] 
energy = [] 
filepath=''
filename=''
with open(filepath+filename, 'r') as f:
    for line in f: # 遍历文件
        if not line: continue 
        t, e = line.split() 
        time.append(float(t)),
        energy.append(float(e))
c00=mean(energy)
time=np.array(time)
popt, pcov = curve_fit(fourier, time, energy)
print(popt)

报错原因与解决方案

报错原因

第一版触发IndexError：curve_fit默认将函数签名中x之后的所有入参作为待拟合的独立标量参数，自动迭代时会将单个数值依次传入。第一版函数将系数数组v、角速度w1、w2设为三个独立入参，curve_fit无法识别v是数组类型，会将第一个标量参数赋值给v，对标量做数组下标索引必然报错。
第二版触发参数数量识别错误：curve_fit要求模型函数第一个入参必须是自变量x，后续入参为待拟合参数。第二版将所有参数（含x）打包为单个args入参，函数签名不符合解析规则，库无法自动统计待拟合参数的数量。
原第二版代码存在笔误：内层循环索引更新写为k*=2，会导致索引跳变，正确的累加逻辑应为k+=2。

修正实现

通过闭包动态生成对应阶数的拟合函数，让生成的函数签名完全匹配curve_fit的要求，同时传入合理的初始参数避免拟合不收敛，可直接运行的修正代码如下：

import numpy as np
from scipy.optimize import curve_fit
from statistics import mean
import matplotlib.pyplot as plt

def build_2d_fourier_model(m, c00):
    """
    动态生成m阶二维傅里叶拟合模型
    返回模型符合curve_fit签名要求：第一个参数为自变量x，后续全为待拟合标量参数
    参数顺序：w1, w2, 后续按顺序排列所有余弦、正弦项系数
    """
    # 预计算所有频率项的索引，避免拟合过程重复计算
    term_indices = []
    # 加入j=0的项，i取值1到m
    for i in range(1, m+1):
        term_indices.append( (i, 0) )
    # 加入j从1到m的项，i取值-m到m
    for j in range(1, m+1):
        for i in range(-m, m+1):
            term_indices.append( (i, j) )
    total_coeff = len(term_indices) * 2 # 每个频率对对应cos、sin两个系数
    total_params = 2 + total_coeff # 总参数=2个角速度+所有级数系数

    def fourier_model(x, *params):
        w1, w2 = params[0], params[1]
        coeffs = params[2:]
        f = np.full_like(x, c00, dtype=np.float64)
        coeff_idx = 0
        for (i,j) in term_indices:
            angle = i * w1 * x + j * w2 * x
            f += coeffs[coeff_idx] * np.cos(angle)
            f += coeffs[coeff_idx+1] * np.sin(angle)
            coeff_idx += 2
        return f
    
    return fourier_model, total_params

# 调用示例
if __name__ == "__main__":
    m = 2 # 可自定义傅里叶级数阶数
    # 替换为实际的数据读取逻辑
    # time, energy = 从txt读取的时序、能量数据
    # 以下为模拟测试数据
    time = np.linspace(0, 10, 1000)
    true_w1, true_w2 = 1.2, 0.8
    true_c00 = 50
    energy = true_c00 + 2*np.cos(1*true_w1*time) + 1.5*np.sin(1*true_w1*time)
    energy += 1*np.cos(-1*true_w1*time + 1*true_w2*time) + 0.8*np.sin(-1*true_w1*time +1*true_w2*time)
    energy += np.random.normal(0, 0.1, size=len(time)) # 加入噪声模拟真实数据

    c00 = mean(energy)
    model, n_params = build_2d_fourier_model(m, c00)
    # 初始化参数：w1、w2建议先通过FFT取频谱峰值作为初始值，系数初始值可设为1
    p0 = np.ones(n_params)
    p0[0], p0[1] = 1.0, 1.0
    popt, pcov = curve_fit(model, time, energy, p0=p0, maxfev=100000)
    print("拟合w1: %.4f, 拟合w2: %.4f" % (popt[0], popt[1]))
    print("拟合系数数量：%d" % (len(popt)-2))

    # 拟合结果可视化
    plt.plot(time, energy, label='原始数据', alpha=0.5)
    plt.plot(time, model(time, *popt), label='拟合曲线', color='red')
    plt.legend()
    plt.show()

注意事项

角速度w1、w2的初始值对拟合结果影响极大，如果初始值离真实值过远，很容易陷入局部最优或拟合不收敛，建议先对输入数据做快速傅里叶变换，取两个最强的频率分量作为初始值。
阶数m不要设置过大，否则参数数量随m呈平方级增长，不仅会大幅提升拟合耗时，还容易出现过拟合问题。

内容的提问来源于stack exchange，提问作者Gilberto

Scipy curve_fit可变参数拟合m阶二维傅里叶级数报错解决

第一版函数报错 #

第二版函数报错 #

报错原因 #

修正实现 #

注意事项 #

第一版函数报错

第二版函数报错

报错原因

修正实现

注意事项