Create scipy curve fitting definitions for fourier series dynamically

thewaywewalk

I'd like to achieve a fourier series development for a x-y-dataset using numpy and scipy.

At first I want to fit my data with the first 8 cosines and plot additionally only the first harmonic. So I wrote the following two function defintions:

# fourier series defintions
tau = 0.045
def fourier8(x, a1, a2, a3, a4, a5, a6, a7, a8):
    return a1 * np.cos(1 * np.pi / tau * x) + \
           a2 * np.cos(2 * np.pi / tau * x) + \
           a3 * np.cos(3 * np.pi / tau * x) + \
           a4 * np.cos(4 * np.pi / tau * x) + \
           a5 * np.cos(5 * np.pi / tau * x) + \
           a6 * np.cos(6 * np.pi / tau * x) + \
           a7 * np.cos(7 * np.pi / tau * x) + \
           a8 * np.cos(8 * np.pi / tau * x)


def fourier1(x, a1):
    return a1 * np.cos(1 * np.pi / tau * x)

Then I use them to fit my data:

# import and filename
filename = 'data.txt'
import numpy as np
import matplotlib.pyplot as plt
from scipy.optimize import curve_fit

z, Ua = np.loadtxt(filename,delimiter=',', unpack=True)
tau = 0.045

# plot data
fig = plt.figure()
ax1 = fig.add_subplot(111)
p1, = plt.plot(z,Ua)

# fits
popt, pcov = curve_fit(fourier8, z, Ua)

# further plots
Ua_fit8 = fourier8(z,*popt)
Ua_fit1 = fourier1(z,popt[0])
p2, = plt.plot(z,Ua_fit8)
p3, = plt.plot(z,Ua_fit1)

plt.show()

which works as desired:

enter image description here

But know I got stuck making it generic for arbitary orders of harmonics, e.g. I want to fit my data with the first fifteen harmonics and plot the first three harmonics.

How could I achieve that without defining fourier1, fourier2, fourier3 ... , fourier15?


Here is the example text file data.txt to play with it:

1.000000000000000021e-03,2.794682735905079767e+02
2.000000000000000042e-03,2.792294526290349950e+02
2.999999999999999629e-03,2.779794770690260179e+02
4.000000000000000083e-03,2.757183469104809888e+02
5.000000000000000104e-03,2.734572167519349932e+02
5.999999999999999258e-03,2.711960865933900209e+02
7.000000000000000146e-03,2.689349564348440254e+02
8.000000000000000167e-03,2.714324329982829909e+02
8.999999999999999320e-03,2.739299095617229796e+02
1.000000000000000021e-02,2.764273861251620019e+02
1.100000000000000110e-02,2.789248626886010243e+02
1.199999999999999852e-02,2.799669443683339978e+02
1.299999999999999940e-02,2.795536311643609793e+02
1.400000000000000029e-02,2.791403179603880176e+02
1.499999999999999944e-02,2.787270047564149991e+02
1.600000000000000033e-02,2.783136915524419805e+02
1.699999999999999775e-02,2.673604939531239779e+02
1.799999999999999864e-02,2.564072963538059753e+02
1.899999999999999953e-02,2.454540987544889958e+02
2.000000000000000042e-02,2.345009011551709932e+02
2.099999999999999784e-02,1.781413355804160119e+02
2.200000000000000219e-02,7.637540203022649621e+01
2.299999999999999961e-02,-2.539053151996269975e+01
2.399999999999999703e-02,-1.271564650701519952e+02
2.500000000000000139e-02,-2.289223986203409993e+02
2.599999999999999881e-02,-2.399383538664330047e+02
2.700000000000000316e-02,-2.509543091125239869e+02
2.800000000000000058e-02,-2.619702643586149975e+02
2.899999999999999800e-02,-2.729862196047059797e+02
2.999999999999999889e-02,-2.786861050144170235e+02
3.099999999999999978e-02,-2.790699205877460258e+02
3.200000000000000067e-02,-2.794537361610759945e+02
3.300000000000000155e-02,-2.798375517344049968e+02
3.399999999999999550e-02,-2.802213673077350222e+02
3.500000000000000333e-02,-2.776516459805940258e+02
3.599999999999999728e-02,-2.750819246534539957e+02
3.700000000000000511e-02,-2.725122033263129993e+02
3.799999999999999906e-02,-2.699424819991720028e+02
3.899999999999999994e-02,-2.698567311502329744e+02
4.000000000000000083e-02,-2.722549507794930150e+02
4.100000000000000172e-02,-2.746531704087540220e+02
4.199999999999999567e-02,-2.770513900380149721e+02
4.299999999999999656e-02,-2.794496096672759791e+02
4.400000000000000439e-02,-2.800761105821779893e+02
4.499999999999999833e-02,-2.800761105821779893e+02
4.599999999999999922e-02,-2.800761105821779893e+02
4.700000000000000011e-02,-2.800761105821779893e+02
4.799999999999999406e-02,-2.788333722531979788e+02
4.900000000000000189e-02,-2.763478955952380147e+02
5.000000000000000278e-02,-2.738624189372779938e+02
5.100000000000000366e-02,-2.713769422793179729e+02
5.199999999999999761e-02,-2.688914656213580088e+02
5.299999999999999850e-02,-2.715383673199499981e+02
5.400000000000000633e-02,-2.741852690185419874e+02
5.499999999999999334e-02,-2.768321707171339767e+02
5.600000000000000117e-02,-2.794790724157260229e+02
5.700000000000000205e-02,-2.804776351435970128e+02
5.799999999999999600e-02,-2.798278589007459800e+02
5.899999999999999689e-02,-2.791780826578950041e+02
5.999999999999999778e-02,-2.785283064150449945e+02
6.100000000000000561e-02,-2.778785301721940186e+02
6.199999999999999956e-02,-2.670252067497989970e+02
6.300000000000000044e-02,-2.561718833274049985e+02
6.400000000000000133e-02,-2.453185599050100052e+02
6.500000000000000222e-02,-2.344652364826150119e+02
6.600000000000000311e-02,-1.780224826854309867e+02
6.700000000000000400e-02,-7.599029851345700592e+01
6.799999999999999101e-02,2.604188565851649884e+01
6.900000000000000577e-02,1.280740698304900036e+02
7.000000000000000666e-02,2.301062540024639986e+02
7.100000000000000755e-02,2.404921248105050040e+02
7.199999999999999456e-02,2.508779956185460094e+02
7.299999999999999545e-02,2.612638664265870148e+02
7.400000000000001021e-02,2.716497372346279917e+02
7.499999999999999722e-02,2.773051723900500178e+02
7.599999999999999811e-02,2.782301718928520131e+02
7.699999999999999900e-02,2.791551713956549747e+02
7.799999999999999989e-02,2.800801708984579932e+02
7.900000000000001465e-02,2.810051704012610116e+02
8.000000000000000167e-02,2.785107135689390248e+02
8.099999999999998868e-02,2.760162567366169810e+02
8.200000000000000344e-02,2.735217999042949941e+02
8.300000000000000433e-02,2.710273430719730072e+02
8.399999999999999134e-02,2.706544464035359852e+02
8.500000000000000611e-02,2.724031098989830184e+02
8.599999999999999312e-02,2.741517733944299948e+02
8.699999999999999400e-02,2.759004368898779944e+02
8.800000000000000877e-02,2.776491003853250277e+02
8.899999999999999578e-02,2.783445666445250026e+02
8.999999999999999667e-02,2.790400329037249776e+02
josliber

I would approach this by defining a single function fourier that accepts a variable number of arguments, using a loop through the cosines to evaluate the function for each of your input points:

def fourier(x, *a):
    ret = a[0] * np.cos(np.pi / tau * x)
    for deg in range(1, len(a)):
        ret += a[deg] * np.cos((deg+1) * np.pi / tau * x)
    return ret

Because this function takes a variable number of arguments, you need to provide a starting position of the correct dimensionality as a hint to the curve_fit function so it knows the number of cosines to use. In the example you provide with 8 cosines:

popt, pcov = curve_fit(fourier, z, Ua, [1.0] * 8)

Now, the use case you're asking about (fitting with 15 harmonics and plotting the first three) can be accomplished with:

# Fit with 15 harmonics
popt, pcov = curve_fit(fourier, z, Ua, [1.0] * 15)

# Plot data, 15 harmonics, and first 3 harmonics
fig = plt.figure()
ax1 = fig.add_subplot(111)
p1, = plt.plot(z,Ua)
p2, = plt.plot(z, fourier(z, *popt))
p3, = plt.plot(z, fourier(z, popt[0], popt[1], popt[2]))
plt.show()

enter image description here

Collected from the Internet

Please contact [email protected] to delete if infringement.

edited at
0

Comments

0 comments
Login to comment

Related

From Dev

Exponential curve fitting in SciPy

From Dev

Scipy sigmoid curve fitting

From Dev

Calculating Fourier series in SciPy

From Dev

python numpy/scipy curve fitting

From Dev

Curve fitting in Python using scipy

From Dev

scipy curve fitting negative value

From Dev

How to calculate the likelihood of curve-fitting in scipy?

From Dev

scipy splrep() with weights not fitting the given curve

From Dev

curve fitting optimisation scipy python numpy

From Dev

Fitting a polinomial curve to time series data

From Dev

Exponential Fitting with Scipy.Optimise Curve_fit not working

From Dev

Correct fitting with scipy curve_fit including errors in x?

From Dev

How can I define the period before fitting a Fourier series to discrete data using MATLAB?

From Dev

Fitting a 2D Gaussian function using scipy.optimize.curve_fit - ValueError and minpack.error

From Dev

Fitting exponential function through two data points with scipy curve_fit

From Dev

Curve Fitting to Exponential

From Dev

Simple curve fitting

From Dev

weighted curve fitting with lsqcurvefit

From Dev

Bayesian curve fitting model

From Dev

Iteratively fitting polynomial curve

From Dev

nonlinear curve fitting in R

From Dev

clustering algorithm with curve fitting

From Dev

Fitting a gaussian to a curve in Python

From Dev

curve fitting with lmfit python

From Dev

Curve Fitting - DataSet

From Dev

Predictive curve fitting matlab

From Dev

Fitting a curve in the points

From Dev

Fitting distribution on curve

From Dev

fitting a distribution to survival curve