DFT

표본 $N$ 개에 어떤 주파수가 얼마나 들어 있는지 어떻게 계산할까?

DTFT와 DFT

표본 수열의 스펙트럼은 연속 신호의 푸리에 변환에 표본화를 이어 붙인 결과이고, 이것이 DTFT $X(\Omega)=\sum_n x[n]e^{-j\Omega n}$ 이다. $\Omega=2\pi f/f_s$ (rad/sample)는 표본 하나마다 도는 각이고, $X(\Omega)$ 는 $2\pi$ 마다 되풀이된다(신호 및 시스템 7장 이산시간 푸리에 변환).

컴퓨터가 가진 것은 표본 $N$ 개이다. 그 $N$ 개와 회전 $N$ 개의 내적이 DFT이다.

X[k]=\sum_{n=0}^{N-1}x[n]\,e^{-j2\pi kn/N}

$n$ 은 시간 지표, $k$ 는 주파수 지표이고 둘 다 $0,1,\ldots,N-1$ 이며, $N$ 은 표본 수이다. 칸 $k$ 는 물리 주파수 $f_k=kf_s/N$ 에 해당하고, 칸 간격은 $\Delta f=f_s/N$ 이다. 관측 시간 $T_{\mathrm{obs}}=N/f_s$ 의 역수와 같다(신호 및 시스템 7장 DFT: 회전과의 내적).

DFT 행렬

$k=0,\ldots,N-1$ 의 DFT를 한꺼번에 쓰면 행렬 곱 $\mathbf X=W\mathbf x$ 이다. $W[k,n]=e^{-j2\pi kn/N}$ 이고, 원소는 모두 1의 $N$ 제곱근이다. 행 $k$ 는 열을 하나 넘길 때마다 $-2\pi k/N$ 씩 도는 시계이다.

아래 그림은 $N=8$ 일 때 $\mathbf X=W\mathbf x$ 를 시계로 그린 것이다. 입력은 $\mathbf x=[1,2,0,-1,0,0,0,0]^{\mathsf T}$ 이다. 먼저 '▶'로 행을 하나씩 내려가며 오른쪽 $X$ 열이 채워지는 것을 보고, 행 2와 행 6의 값을 비교한다.

N = 8 시계 행렬로 쓴 X = Wx이다. 행 k의 시계를 위 줄의 x[n]과 자리마다 곱해 더한 값이 오른쪽 열의 X[k]이다. 신호 및 시스템 7장의 그림을 그대로 썼다.

해 보기

행 2까지 내려간다.
관찰
행 2의 바늘은 열마다 −90°씩 돈다. $X[2]=1+2(-j)-j=1-j3$ 이다.
행 6까지 내려간다.
관찰
$X[6]=1+j3$ 으로 $X[2]$ 의 켤레이다. 행 6의 바늘은 행 2와 같은 속도로 반대 방향으로 돈다.

numpy의 np.fft.fft는 이 행렬 곱과 같은 값을 준다. $N$ 으로 나누지 않으므로, 진폭 1인 복소 톤이 칸에 정확히 놓이면 $|X[k]|=N$ 이다.

N = 8
n = np.arange(N)
W = np.exp(-2j * np.pi * np.outer(n, n) / N)   # W[k, n]
x = np.array([1, 2, 0, -1, 0, 0, 0, 0])
np.allclose(W @ x, np.fft.fft(x))              # True

음의 주파수와 대칭

행 $N-k$ 는 행 $k$ 의 켤레이다. 같은 속도로 반대 방향으로 도는 회전이므로, 칸 $k>N/2$ 는 음의 주파수 $(k-N)f_s/N$ 로 읽는다. $N$ 이 짝수이면 칸 $N/2$ 는 $+f_s/2$ 와 $-f_s/2$ 가 같은 칸이고, numpy는 이 칸을 $-f_s/2$ 로 적는다. 그래서 np.fft.fftshift로 늘어놓은 축은 $-f_s/2\le f<f_s/2$ 이고 $+f_s/2$ 칸은 따로 없다.

예를 들어 $f_s=2$ MHz, $N=400$ 이면 $\Delta f=5$ kHz이다. 칸 60은 $+300$ kHz이고, 칸 340은 $(340-400)\cdot5=-300$ kHz이다.

입력 $x[n]$ 이 실수이면 두 칸의 값은 켤레이다.

\begin{gathered} X[N-k]=\overline{X[k]}\\ |X[N-k]|=|X[k]| \end{gathered}

$X[0]$ 과 ( $N$ 이 짝수일 때) $X[N/2]$ 는 실수이다. 이 대칭은 $\overline{x[n]}=x[n]$ 을 써서 얻으므로, 입력이 복소수이면 성립한다는 보장이 없다. 복소 표본의 스펙트럼이 0 Hz에 대해 대칭이 아닌 이유이다(신호 및 시스템 7장 실수 신호와 복소 신호).

실습에서 쓰는 DFT의 성질

성질 표 전체는 신호 및 시스템 7장 DTFT의 성질에 있다. SDR 코드에서 자주 쓰는 것은 셋이다.

주파수 이동: 표본에 $e^{j2\pi f_0n/f_s}$ 를 곱하면 스펙트럼 전체가 $f_0$ 만큼 옮겨 간다. $f_0$ 이 $\Delta f$ 의 정수배이면 칸으로 $f_0/\Delta f$ 칸 옮겨 가고, 끝을 넘은 칸은 반대쪽 끝으로 돌아온다. 정수배가 아니면 선이 칸 사이에 놓여 이웃 칸으로 번진다(3장). 받은 신호를 소프트웨어로 0 Hz 둘레에 다시 맞추는 방법이 이것이다.
시간 이동: 입력을 $n_0$ 표본 늦추면 $X[k]$ 에 $e^{-j2\pi kn_0/N}$ 이 곱해진다. 크기는 그대로이고 위상만 $k$ 에 비례해 기운다. DFT에서는 끝에서 밀려난 표본이 앞으로 돌아오는 원형 이동이다.
Parseval 관계: $\sum_n|x[n]|^2=\frac1N\sum_k|X[k]|^2$ . 표본 제곱의 합(에너지)을 스펙트럼에서도 같은 값으로 계산한다.