0

Before I get judged, I am no expert with this and simply for the sake of curiosity I tried to write some code that performs a Fourier transform. After watching 3Blue1Brown's video on fourier transform I wanted to write the algorithm myself and plot it, simply because... well it looks cool. I tried to do everything in pure python only using numpy and matplotlib, and it sort of works.

Notice: I am plotting everything iteratively, and re-plotting at each increment

But, plotting the wound up wave and the transform is very very slow. I think I am doing some things inefficiently, maybe even wrong.

Here's what it looks like:

enter image description here

Here is the code:

import numpy as np
import matplotlib.pyplot as plt

u = np.linspace(-8*np.pi, 8*np.pi, 1000)
sin1 = np.sin(u) + 2

u2 = np.linspace(-3*np.pi, 3*np.pi, 1000)
sin2 = np.sin(u2) + 2
plt.plot(u + u2, sin1 + sin2)
fig, (winder, integax) = plt.subplots(nrows = 2, ncols = 1)

L = len(sin1)
real = []
imag = []
integral = []

for val in np.arange(0.00001,360,0.00001):
    real = []
    imag = []
    for t,si in zip(np.arange(0,L,val),sin1 + sin2):
        complex = si * np.e ** (2 * np.pi * 1j * t)
        real.append(complex.real)
        imag.append(complex.imag)

    fig.set_size_inches(10,10)
    point = np.trapz(real)
    integral.append(point)
    #print(integral[-1], time[-1])
    integax.plot(integral)
    winder.plot(real, imag, 'b-')
    plt.pause(0.00001)
    winder.cla()

ax = plt.plot(real, imag, 'b-')
plt.show()

Now I would like to plot it faster, and I think that the integration part is not correct. Since no spikes occur in the resulting plot even after waiting for a long time.

I also don't think that I am using linspace correctly to plot sine waves nor am I doing the frequency part right in the fourier formula.

Ahmad Moussa
  • 876
  • 10
  • 31
  • 1
    If I understood your code, you are using numerical integration with the trapezoidal rule. If you want to compute a Fourier Transform and do it fast, then you should probably implement a [Fast Fourier transform](https://en.wikipedia.org/wiki/Fast_Fourier_transform). Note that the FFT is considered (along the Kalman Filter) to be one of the most important algorithms in engineering. – Mefitico Jun 10 '19 at 12:08
  • I have heard of the fast fourier transform and wanted to tackle it after I mastered the discrete, as I'm afraid that it might be a little be too advanced for me. In any case, I can live with the speed constraints right now, as this is just a little demonstration I wanted to present for my lab. I just that maybe the slowness could be due to how I plot the wave and the integral of it, and not the integration itself. – Ahmad Moussa Jun 10 '19 at 13:57
  • 1
    You're approximating a continuous Fourier transform - the discrete Fourier transform(DFT) would use a sum rather than integration. – user2699 Jun 10 '19 at 14:04
  • 2
    @AhmadMoussa : If you are not sure about which step is taking longer in your code, then you need to perform a profiling operation on it, see for instance [this question](https://stackoverflow.com/questions/582336/how-can-you-profile-a-python-script) – Mefitico Jun 10 '19 at 14:09
  • 2
    `matplotlib` can introduce quite a bit of overhead when you repeatedly plot (although @Metifico is right about profiling to be sure). Take a look at https://stackoverflow.com/questions/48697184/faster-plotting-of-real-time-audio-signal/48698434#48698434 for some hints on how to speed up matplotlib (although it's about `imshow`, `plot` can be treated similarly). – user2699 Jun 10 '19 at 14:11
  • Thank you both for the suggestions, I will make sure to have a look at the link. – Ahmad Moussa Jun 10 '19 at 14:20

1 Answers1

0

Compiling comments into an answer:

If I understood your code, you are using numerical integration with the trapezoidal rule. If you want to compute a Fourier Transform and do it fast, then you should probably implement a Fast Fourier transform algorithm. Note that the FFT is considered (along with the Kalman Filter) to be one of the most important algorithms in engineering, and for a good reason: It is so much faster than naive implementations of the Fourier transform that it turns intractable problems into computing demanding but feasible products.

As mentioned be @user2699, there are methods to speed up plotting, as per this question.

This addresses some algorithmic and Python related bottlenecks of your code, but if you are dissatisfied with the speed it takes to run, the go-to starting point should be profiling your code, and some hints on how to do this with Python have been discussed in this question.

Mefitico
  • 816
  • 1
  • 12
  • 41