Fourier Series

Try to calculate a fourier series in python



In [1]:

    
from __future__ import division
import numpy as np
from scipy import integrate
import matplotlib.pyplot as plt

main functions

use sn = fourier_series(f, N, offset=0, period=2*np.pi) to create a fourier series of function $f(x)$ in range $(offset,\; offset+period)$
use sn(t) to calculate $S_n(x)$



In [2]:

    
class fourier_series :
    def __init__(self, f, N, offset=0, period=2*np.pi) :
        self.N = N
        self.period = period
        self.offset = offset
        self.c = np.zeros(2*N+1, dtype=complex)
        for n in range(-N, N+1) :
            self.c[n] = complex(integrate.quad(lambda x : (f(x)*np.exp(-2.0j*np.pi*n*x/period)).real, \
                                          offset, offset + period)[0], \
                           integrate.quad(lambda x : (f(x)*np.exp(-2.0j*np.pi*n*x/period)).imag, \
                                          offset, offset + period)[0])/period
        self._v_sn = np.vectorize(self._sn)
    
    def _sn(self, x) :
        ans = 0.0j
        for n in range(-self.N, self.N+1) :
            ans += self.c[n] * np.exp(2.0j*np.pi*n*(x-self.offset)/self.period)
        return ans
    
    def __call__(self, x) :
        return self._v_sn(x)

Example 1:



In [3]:

    
@np.vectorize
def f(x) :
    return np.floor(x)%2

s = fourier_series(f, 10, 0., 2.)
s









    Out[3]:





<__main__.fourier_series instance at 0x11247cfc8>



In [4]:

    
t = np.linspace(0, 4, 200)
plt.plot(t, f(t), 'rx', t, s(t).real, 'b-')
plt.legend(['original', 'fourier'])
plt.title('Fourier Series --Example 1')
plt.show()

Another example

Here we want to calculate the function which is given by a series of $\sin$s and $\cos$s of a star.



In [5]:

    
from scipy.interpolate import interp1d

To get started, at first we find the vertices in complex value like this:



In [6]:

    
s = np.array([np.exp(4.0j*np.pi*x/5.) for x in range(6)])
s









    Out[6]:





array([ 1.00000000 +0.00000000e+00j, -0.80901699 +5.87785252e-01j,
        0.30901699 -9.51056516e-01j,  0.30901699 +9.51056516e-01j,
       -0.80901699 -5.87785252e-01j,  1.00000000 -4.89858720e-16j])

Then, we may transfer the lines into a $\mathbb{R} \rightarrow \mathbb{C}$ function:



In [7]:

    
_t = np.linspace(0., 2*np.pi, num=len(s))
g = interp1d(_t, s)

And now create a Fourier Series:



In [8]:

    
fs = fourier_series(g, 4)
t = np.linspace(0., 2*np.pi, 200)
fs.c









    Out[8]:





array([ -7.74430137e-17 -4.20238982e-18j,
         8.77521433e-17 +8.71966840e-17j,
         5.72786697e-01 -7.34546505e-17j,
        -1.61488402e-14 +2.96562390e-17j,
         5.45145663e-14 -1.02242009e-16j,
        -5.96152716e-14 +2.38382337e-16j,
         2.54571865e-01 -2.73638464e-16j,
         1.60049634e-15 -1.51945706e-16j,   4.09369540e-15 -1.12804325e-17j])



In [9]:

    
plt.plot(fs(t).real, fs(t).imag, 'b-', label='fourier')
plt.plot(g(t).real, g(t).imag, 'g--', label='interpolate')
plt.plot(s.real, s.imag, 'bx', label='original')
plt.legend()
plt.title('Fourier Series --star')
plt.show()