Блог

Я пытаюсь найти ряд Фурье, который соответствует этому графику:

Сначала я воспользовался этим сайтом и нашел следующие коэффициенты Фурье

 a = [0.983, 0.292, 0.110, 0.006, -0.022, -0.018, 0.006]
b = [0.000, 0.246, 0.002, -0.021, 0.017, 0.046, 0.019]
omega = 0.045 
  

что привело к следующему результату

Это довольно хорошо, но затем я решил попробовать закодировать свой собственный curve fitter на python и обнаружил, что библиотека symfit кажется лучшим вариантом. Я даже переделал один пример, найденный в документации, изменив очень мало.

И, конечно же, когда я запускаю код, подгонка кажется довольно хорошей на графиках, созданных непосредственно с помощью объекта fit:

plt.plot(xdata, fit.model(x=xdata, **fit_result.params).y, color='green', ls=':')

Но затем я решил использовать коэффициенты Фурье напрямую, а не полагаться на объект fit, чтобы попытаться воспроизвести функцию. Вот как код печатает эти коэффициенты:

 Parameter Value        Standard Deviation
a0        -1.206104e 03 nan
a1        1.799223e 03 1.087713e 05
a10       2.393892e 00 8.023119e 02
a2        -8.465307e 02 2.133670e 05
a3        6.798590e 02 nan
a4        -9.287399e 02 nan
a5        7.388156e 02 nan
a6        -2.099387e 02 1.376980e 05
a7        -9.838123e 01 9.799211e 04
a8        9.821148e 01 3.739932e 04
a9        -2.827002e 01 8.152094e 03
b1        -6.084419e 02 2.480938e 05
b10       4.586827e 01 3.396681e 02
b2        9.311030e 02 2.178771e 05
b3        -8.877734e 02 2.406878e 05
b4        8.810145e 02 4.462669e 05
b5        -1.288586e 03 4.376938e 05
b6        1.734395e 03 2.730319e 05
b7        -1.583775e 03 1.137506e 05
b8        9.115344e 02 3.084382e 04
b9        -3.045882e 02 4.905459e 03
w         -1.370387e-02 nan
  

Стандартные отклонения были огромным красным флагом, они абсурдно огромны, и, конечно же, это результат

Это похоже на плохую шутку.

Это было сделано с использованием 10 коэффициентов, если я использую только шесть, я получаю что-то едва приемлемое.

Хуже всего то, что графики, созданные непосредственно с помощью объекта fit, должны означать, что код действительно хорошо подходит, но по какой-то причине коэффициенты Фурье, которые он мне дает, просто ложь.

Вот мой код.

Во-первых, вот как я подстраиваю функции

 from symfit import parameters, variables, sin, cos, Fit
import numpy as np
import matplotlib.pyplot as plt
import math as m

def fourier_series(x, f, n=0):
    a0,*cos_a = parameters(','.join(['a{}'.format(i) for i in range(0, n   1)]))
    sin_b = parameters(','.join(['b{}'.format(i) for i in range(1, n   1)]))
    series = a0   sum(ai * cos(i * f * x)   bi * sin(i * f * x)
                     for i, (ai, bi) in enumerate(zip(cos_a, sin_b), start=1))
    return series

x, y = variables('x, y')
w, = parameters('w')
model_dict = {y: fourier_series(x, f=w, n=10)}
print(model_dict)

#This is the graph I'm trying to fit turned into data points
xlist=[381.854,385.274,387.55,388.973,390.656,392.671,394.281,396.37,398.014,402.047,408.288,411.885,414.452,418.099,
       423.493,426.901,428.014,430.068,433.114,438.699,445.023,448.562,453.307,458.014,463.663,467.877,474.536,
       479.384,485.409,489.247,496.283,501.164,508.191,514.315,520.1]
ylist=[0.02816,0.09827,0.18838,0.2948,0.36316,0.46243,0.53795,0.60116,0.64162,0.68361,0.71676,0.79285,0.87861,0.94579,
       0.99422,0.95307,0.89017,0.85549,0.79285,0.76301,0.7273,0.65896,0.59621,0.53757,0.47969,0.43353,0.38501,
       0.33526,0.28305,0.23121,0.18838,0.14451,0.11555,0.08092,0.07185]

xdata = np.array(xlist)
ydata = np.array(ylist)

# Define a Fit object for this model and data
fit = Fit(model_dict, x=xdata, y=ydata)
fit_result = fit.execute()
print(fit_result)

# Plot the result
plt.plot(xdata, ydata)
plt.plot(xdata, fit.model(x=xdata, **fit_result.params).y, color='green', ls=':')
  

Вот как я беру эти коэффициенты и пытаюсь воспроизвести график. У меня есть разные наборы данных a b и omega, чтобы попытаться воспроизвести другие ряды Фурье. Единственный, который работает, — это шаговый клин (пример в документации).

 from symfit import parameters, variables, sin, cos, Fit
import numpy as np
import matplotlib.pyplot as plt
import math as m
from random import random

def f(x,a,b,omega):
    output = a[0]/2
    for k in range(1,len(a)-1):
        output = output   a[k]*m.cos(k*omega*x)   b[k]*m.sin(k*omega*x)
    return output

def monteCarlo(min,max,a,b,omega):
    x = min   random()*(max-min)
    return (x,f(x,a,b,omega))

xdata=[]
ydata=[]

#for step function between -3 and 3
#a=[5.000000e-01,1.106075e-10,2.364473e-11,1.697959e-10,2.822410e-11,3.860925e-11,-1.522386e-10]
#b=[0,6.267589e-01,2.020945e-02,1.846406e-01,3.623074e-02,8.726419e-02,4.518721e-02]
#omega=8.615115e-01

#for our spectrum with n=6
a = [-1.755069e 02,-1.097051e 02,1.571135e 02,6.491389e 01,1.220463e 02,1.778165e 02,4.564577e 01]
b = [0,-5.523502e 01,1.517397e 02,1.325198e 02,-1.015198e 02,-3.255781e 01,2.495038e 01]
omega = -5.353751e-03

#for our spectrum with n=10
#a = [-1.206104e 03,1.799223e 03,-8.465307e 02,6.798590e 02,-9.287399e 02,7.388156e 02,-2.099387e 02,-9.838123e 01,
#    9.821148e 01,-2.827002e 01,2.393892e 00]
#b = [0,-6.084419e 02,9.311030e 02,-8.877734e 02,8.810145e 02,-1.288586e 03,1.734395e 03,-1.583775e 03,9.115344e 02,
#     -3.045882e 02,4.586827e 01]
#omega = -1.370387e-02

#a = [0.983, 0.292, 0.110, 0.006, -0.022, -0.018, 0.006]
#b = [0.000, 0.246, 0.002, -0.021, 0.017, 0.046, 0.019]
#omega = 0.045

for n in range(0,1000):
    coor = monteCarlo(375,525,a,b,omega)
    xdata.append(coor[0])
    ydata.append(coor[1])

plt.plot(xdata, ydata, 'bo')