Matplotlib로 FFT 그려보기

Matplotlib를 사용하여 FFT를 그려보았습니다. 

2021. 12. 18  최초작성

2021. 12. 31  예제 코드 추가

2021. 01. 17  예제 코드 추가 

아직 정확히 개념이 잡힌 상태가 아니라서 틀린점이 있을 수 있으니 참고용으로만 사용하세요.

다음 두 개의 sin 그래프를 FFT로 변환하여 그려봅니다.

오른쪽의 진폭이 2배 더 큽니다.  

왼쪽은 FFT로 변환한 결과이며 오른쪽은 FFT 결과를 역변환하여 얻은 원래 파형을 그린것입니다.

두 개의 sin 그래프의 FFT 결과에서 다른 점은 1Hz에서의 높이가 원래 파형의 배수만큼 다른 것입니다.

원래 그래프에서 2배 차이가 났었는데 FFT의 1Hz에서의 높이도 2배 차이입니다.  

 

 

전체 소스코드입니다.

 

# 참고 : https://pythonnumericalmethods.berkeley.edu/notebooks/chapter24.04-FFT-in-Python.html import numpy as np from matplotlib import pyplot as plt from numpy.fft import fft, ifft x = np.arange(0, 60) y1 = np.sin(7*x * np.pi / 180. ) y2 = np.sin(7*x * np.pi / 180. + 3 ) * 2 plt.figure(figsize = (12, 6)) plt.subplot(121) plt.plot(x, y1, color='r') plt.subplot(122) plt.plot(x, y2, color='b') plt.show() data = [y1, y2] for d in data:     x = d     sr = 60 # sampling rate     ts = 1.0/sr # sampling interval     t = np.arange(0, 1, ts)     X = fft(x)     N = len(X)     n = np.arange(N)     T = N/sr     freq = n/T     plt.figure(figsize = (12, 6))     plt.subplot(121)     plt.stem(freq, np.abs(X), 'b', markerfmt=" ", basefmt="-b")     plt.xlabel('Freq (Hz)')     plt.ylabel('FFT Amplitude |X(freq)|')     plt.xlim(0, 10)     plt.subplot(122)     plt.plot(t, ifft(X), 'r')     plt.xlabel('Time (s)')     plt.ylabel('Amplitude')     plt.tight_layout()     plt.show()

sine 파를 생성할 때 주파수로 7 Hz를 사용하였더니 FFT 결과 7 HZ 위치에만 값이 존재합니다. 

 

전체 코드입니다.

 

# 참고 # https://tuttozurich.tistory.com/15 # https://github.com/hazevt04/fft_sine_wave/blob/main/fft_sine_wave.py import numpy as np from matplotlib import pyplot as plt freq = 7          # Frequency of the Sinewave rate = 1000       # Sampling Rate time_interval = 3 # Sampling Interval t = np.linspace(0, time_interval, time_interval*rate) # Time-Axis y1 = 2*np.sin(2*np.pi*freq*t) # Sinewave y2 = 4*np.sin(2*np.pi*freq*t + 2) # Sinewave def PositiveFFT(x, Fs):     Ts = 1/Fs     Nsamp = x.size     xFreq = np.fft.rfftfreq(Nsamp, Ts)[:-1]     yFFT = (np.fft.rfft(x)/Nsamp)[:-1]*2     return xFreq, yFFT # Calculate FFT xFreq1, FFT_data1 = PositiveFFT(y1, rate) xFreq2, FFT_data2 = PositiveFFT(y2, rate) # Plot the FFT fig, ((ax1, ax2), (ax3, ax4)) = plt.subplots(2, 2, figsize = (10,5)) ax1.plot(t, y1) ax1.set_title('signal') ax1.set_xlabel('time') ax1.set_ylabel('amplitude') ax1.set_xlim(0,1) ax2.stem(xFreq1, abs(FFT_data1), use_line_collection=True) ax2.grid(True) ax2.set_title(r"$\frac{2 \cdot fft(s)}{len(s)}$") ax2.set_xlabel('freq') ax2.set_ylabel('amplitude') ax2.set_xlim(0,10) ax3.plot(t, y2) ax3.set_title('signal') ax3.set_xlabel('time') ax3.set_ylabel('amplitude') ax3.set_xlim(0,1) ax4.stem(xFreq2, abs(FFT_data2), use_line_collection=True) ax4.grid(True) ax4.set_title(r"$\frac{2 \cdot fft(s)}{len(s)}$") ax4.set_xlabel('freq') ax4.set_ylabel('amplitude') ax4.set_xlim(0,10) plt.tight_layout() plt.show()

5 Hz, 7 Hz, 9 Hz 사인파를 더한 후, FFT 변환하여 스펙트럼을 그려보았습니다.  사인파를 위한 주파수별로 성분이 검출되었습니다. 

 

# 참고  # https://tuttozurich.tistory.com/15 # https://github.com/hazevt04/fft_sine_wave/blob/main/fft_sine_wave.py # https://alice-secreta.tistory.com/22 # https://techreviewtips.blogspot.com/2017/11/05-02-fft.html import numpy as np from matplotlib import pyplot as plt freq1 = 5   # Frequency of the Sinewave 5 Hz freq2= 7    # 7 Hz freq3 = 9  # 9 Hz Fs = 1000       # Sampling Rate time_interval = 1 # Sampling Interval t = np.arange(0, time_interval, 1/Fs) # 5Hz에 진폭이 2인 시그널을 1초동안 샘플링한다.  y1 = 2*np.cos(2*np.pi*freq1*t) # Sinewave y1 += 2*np.cos(2*np.pi*freq2*t) # Sinewave y1 += 2*np.cos(2*np.pi*freq3*t) # Sinewave def PositiveFFT(x, Fs):     Ts = 1/Fs     Nsamp = x.size     xFreq = np.fft.rfftfreq(Nsamp, Ts)[:-1]     yFFT = (np.fft.rfft(x)/Nsamp)[:-1]*2     #  yFFT = (np.fft.fft(x)/Nsamp)     # yFFT = np.abs(yFFT[:int(Nsamp/2)])*2# 두배를 해야 하는가?     # xFreq = np.fft.fftfreq(Nsamp, Ts)[:int(Nsamp/2)]     return xFreq, yFFT # Calculate FFT xFreq1, FFT_data1 = PositiveFFT(y1, Fs) print(y1.shape, FFT_data1.shape) # Plot the FFT fig, (ax1, ax2) = plt.subplots(1, 2, figsize = (10,5)) ax1.plot(t, y1) ax1.set_title('signal') ax1.set_xlabel('time') ax1.set_ylabel('amplitude') ax1.set_xlim(0,1) ax2.stem(xFreq1, abs(FFT_data1), use_line_collection=True) ax2.grid(True) ax2.set_title(r"$\frac{2 \cdot fft(s)}{len(s)}$") ax2.set_xlabel('freq') ax2.set_ylabel('amplitude') ax2.set_xlim(0,10) plt.tight_layout() plt.show()