파이썬 DFT 코드

미리보기: 결과 파형

이 코드는 Python으로 작성된 이산 푸리에 변환(DFT, discrete Fourier transform) 예제입니다. 주어진 시간축에서 주파수 성분의 크기 및 위상을 추출하기 위해 이산 푸리에 변환을 사용합니다.

코드의 첫 부분에서는 NumPy와 Matplotlib 라이브러리를 import 합니다. DFT 함수는 Numpy 배열을 사용하며, 결과를 그래프로 시각화하는 데에는 Matplotlib가 사용됩니다. (2023.03.23 - [컴퓨터/파이썬] - 파이썬 matplotlib 간단 소개)

import numpy as np
import matplotlib.pyplot as plt

DFT 함수의 입력으로는 N(샘플링 포인트의 수)와 sig(시그널)을 받습니다. 이 함수에서는 N 크기의 윈도우를 정의하고, 이 윈도우를 이용해 입력 시그널 sig를 샘플링합니다. 그리고 각 샘플에서 윈도우를 한 칸씩 이동시키면서 DFT 결과를 계산합니다. 결과는 각 주파수 성분의 크기와 위상값으로 이루어진 배열로 반환됩니다.

def DFT(N, sig):
    window = np.zeros(N)
    component = np.zeros(len(sig), dtype=np.complex128)

    for step in range(len(sig)):
        window[N-1] = sig[step]

        component[step] = (2/N) * np.sum(window * np.exp(-1j*2*np.pi*np.arange(N)/N))

        window[0:N-1] = window[1:N]

    mag = np.abs(component)
    re = np.real(component)
    im = np.imag(component)
    ang = np.arctan2(re, im)
    
    return mag, ang

참고: DFT 수식

다음으로, 주파수와 시간축을 정의합니다. 이 예제에서는 기본주파수 f = 60Hz, N = 128(한 주기당 128개의 샘플), fs = f * N(샘플링 주파수), T = 1/f(주기)로 설정합니다.

# Define the sampling frequency and duration of the signal
N = 128 # sample number per period
f = 60  # fundamental frequency
fs = f * N  # sampling frequency
T = 1/f  # period

그리고 A = 1인 60Hz의 사인파를 생성하고, DFT 함수를 사용하여 이 시그널의 주파수 성분을 계산합니다. 이후 시그널과 DFT 결과를 그래프로 그립니다.

# Generate a sinusoidal wave with frequency f and amplitude A
A = 1
t = t = np.arange(0, 5*T, 1/fs) # time axis
signal = A * np.sin(2*np.pi*f*t)  # sinusoidal signal

# Compute the DFT of the signal using the previously defined function
mag, ang = DFT(N, signal)

# Plot the original signal and its magnitude spectrum
fig, (ax1, ax2) = plt.subplots(nrows=2, ncols=1, figsize=(8, 6), sharex=True)

그래프는 subplot() 함수를 사용하여 2행 1열로 나뉘어 그려집니다. 첫 번째 그래프는 시그널과 주파수 성분의 크기(magnitude)를 보여주며, 두 번째 그래프는 시그널과 주파수 성분의 위상값(angle)을 보여줍니다.

ax1.plot(t, signal, label='Signal')
ax1.plot(t, mag, label='Magnitude')
ax1.set_ylabel('Amplitude')
ax1.set_title('DFT Result - Magnitude')
ax1.legend()
ax1.grid()

ax2.plot(t, signal, label='Signal')
ax2.plot(t, ang, label='Angle')
ax2.set_ylabel('Angle (radians)')
ax2.set_yticks([-np.pi, -np.pi/2, 0, np.pi/2, np.pi])
ax2.set_yticklabels(['$-\pi$', '$-\pi/2$', '$0$', '$\pi/2$', '$\pi$'])
ax2.set_title('DFT Result - Angle')
ax2.legend()
ax2.grid()

plt.tight_layout()
plt.show()

따라서 이 예제는 DFT를 사용하여 시간축에서 주파수 성분을 추출하고, 이를 그래프로 시각화하는 방법을 보여줍니다. 이는 음성인식, 이미지 및 비디오 압축, 신호 처리 및 분석 등 다양한 분야에서 활용됩니다.