• Cute Fourier Transformation coming !

Fourier transformation

$X(t)$ is the amplitude-time function, we want to transform this into frequency domain. Typically, we want to decompose the function into several sine and cosine function with different frequency, phase and amplitude. $F$ represents the frequency we focus on.

$$X(F) = \int_{-\infty}^{+\infty} X(t) e ^{-i2\pi Ft}dt$$

The dot product verifies the similarity of the analysis function and the amplitude-time function.

Discrete Fourier transformation

When we can only sample data from the signals, we replace the X with the folowwing.

$$X_k = \sum_{n=0}^{N-1} X(t) e ^{-\frac{i2\pi kn}{N}}dt$$

$$F=\frac{k}{n}, t=n$$

Expansion of DFT Euler’s equation

$$e^{ix} = cosx + isinx$$

According to Euler’s equation, and regards the $\frac{2 \pi kn}{N}$as $b_n$, we have the following

$$X_k = x0[cos(-b0)+isin(-b0)]+ \dots$$

$$X_k =A_k + B_k i$$

What is A and B w.r.t. k? A and B is the amplitude and phase of cosine wave.

Experiment

Let’s take an experiment!

We first generate a sine signal from 0 to 2 second. We sample the signal with frequency 512 Hz.

# sample frequency is 512 Hz, i.e. sample 512 data per second  
fs = 512  
# time duration 2s  
t = 10  
# number of sample  
N = int(t*fs)  
# sample spacing  
T = 1 / fs  
# generate signals
samples = np.linspace(0,t,N,endpoint=False)  
signals = np.sin(samples) 
  1. Take Fourier Transformation(fft). It returns an array with A/2 +B/2 * i for each frequency from negative to positive.

    yf = fft(signals)
    
  2. Take fftfreq to get the frequency bin. It returns an array with positive frequency bin and negative frequency bin.

    xf = fftfreq(N, T)
    
  3. According to Nyquist limit Nyquist limit = $fs/2$, we only take the positive frequency. We multiply the value in positive frequency by 2 and take fraction of sample number N to get the true FFT value.

    result = 2.0/N * np.abs(yf[0:N//2])
    

The complete code is as follows:

  
from scipy.fftpack import fft, fftfreq  
import numpy as np  
import matplotlib.pyplot as plt  
  
## data generation  
# sample frequency is 512 Hz, i.e. sample 512 data per second  
fs = 512  
# time duration 2s  
t = 10  
# number of sample  
N = int(t*fs)  
# sample spacing  
T = 1 / fs  
# generate signals
samples = np.linspace(0,t,N,endpoint=False)  
signals = np.sin(samples)  
plt.figure()  
plt.plot(samples,signals)  
plt.xlabel("Time(s)")  
plt.ylabel("Amplitude")  
plt.show()  
## fourier transformation  
yf = fft(signals)  
## frequency bin  
"""  
f = [0, 1, ..., n/2-1, -n/2, ..., -1] / (d*n)         if n is even  
f = [0, 1, ..., (n-1)/2, -(n-1)/2, ..., -1] / (d*n)   if n is odd  
"""  
xf = fftfreq(N, T)[:N//2]  
  
# nyquist limit  
nl = fs / 2  

# plot
plt.figure()  
plt.xlabel("Frequency(Hz)")  
plt.ylabel("Amplitude")  
plt.plot(xf,2.0/N * np.abs(yf[0:N//2]))  
plt.show()

We could see the result

great reference

Simon Xu

Wolfram