Fourier series in Python

Viewed 48

I am trying to implement Complex Exponential Fourier Series for f(x) defined on [-L,L] using these formulas,

enter image description here

I want to be able to implement these without calling the Fourier functions in other libraries since I want to also understand what's going on. Here is my attempt,

import numpy as np
from matplotlib import pyplot as plt

steps = 100
dt = 1/steps
L = np.pi
t = np.linspace(-L, L, steps)


def constant(X, Y, n):
    return (1/(2*L))*sum([y*np.exp((1j*n*np.pi*t)/L)*dt for t, y in zip(X, Y)])


def complex_fourier(X, Y, N):
    _X, _Y = [], []
    for t in X:
        f = 0
        for n in range(-N//2, N//2 + 1):
            c = constant(X, Y, n)
            f += c*np.exp((-1j*n*np.pi*t)/L)
        _X += [f.real]
        _Y += [f.imag]
    return _X, _Y


X, Y = complex_fourier(t, np.sin(t), 50)
plt.plot(X, Y, 'k.')
# plt.plot(t, np.sin(t))
plt.show()

The plot seems to be almost random and does not improve with more c terms. Could someone point out exactly what am I doing wrong?

1 Answers

Just to answer the "random plot" part of the question for now - note the Y-scale of your plot!

>>> np.min(Y), np.max(Y)
(-6.1937063114043705e-18, 6.43002899658067e-18)
>>> np.min(X), np.max(X)
(-0.15754356079010426, 0.15754356079010395)

In other words, all of your coefficients are basically real valued. You probably wouldn't be interested in an plot of the imaginary part vs the real part, but rather the sum of squares vs the frequency or mode number.

Related