신호처리에서 다시 만난 FFT: 파형 합성과 하드웨어 구현 고민

예전 이 도메인 eunchan.kim으로 이전하면서 실수로 몇몇 글을 이전하지 못했었습니다. 다른 건 다 괜찮지만, 가장 아쉬운 글이 FFT 관련 글이었죠. 대학교 때에 이해했던 내용을 바탕으로 예제 C 코드까지 있었습니다. 다시 작성해 볼까 싶었지만, 다시 쓰려니 FFT를 싹 잊어먹었더라구요. 그 글을 쓰기 위해 그 글이 필요한 상황이 되었습니다.

이번에 일을 하게 되면서 다시 FFT를 접할 일이 생겼습니다. 이 글에서 신호처리에서 사용되는 FFT에 대해 살펴볼까 합니다.

FFT의 정의 및 동작 원리

FFT를 단순하게 정의하면, 시간 영역의 신호에서 각 주파수의 정보가 얼마나 있는지 찾아서 주파수 영역의 신호로 바꾸는 과정입니다. 반대로 이 주파수 정보를 역변환(iFFT, Inverse FFT)에 인가하면 원래의 시간 영역 신호로 다시 복원할 수 있습니다.

아래의 시뮬레이터는 제가 이전 글에서 그렸던 sine wave, cosine wave를 이용해 어떻게 블랙박스 안에서 만들어지는 신호를 검출하는지를 제미나이에게 구현해달라고 부탁해서 만든 시뮬레이터입니다. 그땐 단순한 그림이었는데, 이젠 검침 주파수를 변경해가면서 어느 주파수에서 잡히는 지를 명확하게 보여줍니다.

처음 제가 FFT를 배울때에는 이 주파수 영역을 2D 이미지에서 이용했었습니다. 이미지의 노이즈를 제거하는 (보통은 노이즈는 고주파 작은 진폭에 몰려있으니까요) 필터를 적용한다던지요.

그 이후로 FFT를 접할 기회가 없었죠. 그러다 최근에 일하면서 신호처리를 다룰 일이 많아지고, 자연스럽게 FFT를 사용할 일이 생겼네요. 이번에는 영상처리가 아니라 파형 신호처리 부분에서 살펴볼까 합니다.

위의 시뮬레이터는 연속적인 주파수 영역을 가정하고 구현한 것입니다. 하지만 실제로 디지털 영억에서는 어떻게든 구간을 쪼개야 하죠. 이미지에서 한 픽셀 한 픽셀 값이 있는것처럼, 신호에서도 한 시각에 하나의 값을 가지고, 다음 시각과 시간이 분리되어 있습니다.

그래서 Discrete Fourier Transform (DFT)가 정의됩니다. 일단 모든 길이를 한번에 변환 하기엔 너무 길고, 시간도 오래 걸리니, 정해진 구간 N 만큼 쪼개서 변환합니다. 만일 신호가 규칙적으로 1, 0, 1, 0 반복한다고 합니다. 이건 N 구간에서 가장 빠르게 회전하는 신호가 될겁니다. 반대로 1, 1, 1, 1 반복한다면 이건 0번 회전하는 신호가 되겠죠. N 구간 안에서 1 -> 0 -> 1 을 cosine wave로 움직인다면, 이건 한번 회전하는 신호가 되죠.

즉, FFT 블럭 사이즈 N 에서는 0번 회전부터 N/2 번 회전까지 표현할 수 있습니다. 회전을 표현하기 위해 FFT에서는 복소평면 (Complex plain)을 사용합니다.

$$ X_k = \sum_{n=0}^{N-1} x_n e^{-j \frac{2\pi k n}{N}} $$

$ e^{-j \frac{2\pi k n}{N}} = cos(\frac{2\pi k n}{N}) - j sin(\frac{2\pi k n}{N}) $ 이죠. 즉, 각 구간에서 현재의 각도에 대한 cosine, sine 값을 이용하면 이 구간에서의 신호에 대한 주파수를 계산할 수 있습니다. 그리고 이 값을 모든 구간에 대해 다 더하면, N 구간에 대한 주파수 영역의 값을 계산할 수 있습니다.

Radix-2

이를 효율적으로 계산하게 하는 방법이 Radix-2 FFT 입니다. Butterfly 계산 (실제로 그림이 나비모양처럼 생겼습니다) 을 단계별로 수행하면 위의 수식을 그대로 계산해 낼 수 있죠. 그래서 Fourier Transform의 복잡도가 $O(N^2)$ 인데, 이를 $O(N log N)$ 으로 줄일 수 있습니다.

FFT Block Size, Frequency Resolution

저는 복잡한 수식을 바로 보면 이해하기 어렵더라구요. 제가 접했던 문제는, 주파수 도메인 신호를 시간 도메인으로 복원하는 부분이었습니다. 이미지 처리만 해봤던 저에게는 그냥 Discrete iFFT를 하면 되는거라고 생각했는데, 이게 딱 명확하게 이해가 안되더라구요.

예를 들어봅시다. 로직이 100MHz로 동작하면 Nyquist 이론에 의해서 최대 50MHz의 신호까지 만들어낼 수 있습니다. 그러면 10MHz 부터 50MHz까지 1MHz 간격의 신호를 만들어내기 위해서는 FFT 블럭사이즈는 몇 개가 필요할까요?

전 처음엔 $\frac{50\text{MHz}}{1\text{MHz}} = 50$ 이니까 50개 block size 라고 생각했습니다. 주파수 성분을 담는 각각의 방인 하나의 bin(빈)이 하나의 주파수를 담당하니까요.

그렇다면 10MHz 부터 50MHz까지 1Hz 간격으로 주파수를 만드려면요? 이땐 block size가 5천만개가 필요할까요? 그렇게 가능은 하겠죠. 하지만, 이걸 소프트웨어나 하드웨어로 구현하려면, 5천만개의 값을 butterfly 곱셈을 해야하니, 긴 시간이 걸리거나, 엄청난 양의 리소스를 사용해야합니다.

제가 바보였더라구요. 항상 이미지의 값을 변환해서 사용하다보니, 딱 떨어지는 값만 상상했나 봅니다. 사실, 주파수와 상관없이 어떤 값이든, 어떤 사이즈의 FFT를 이용해서 신호를 복원할 수 있습니다. 이걸 반대로 특정 주파수의 sine wave가 FFT를 거치면 각 bin에 어떻게 퍼져서 나타나는지를 보면 이해하기 쉬워집니다.

끝없이 이어지는 사인파를 정해진 길이(N개의 샘플)로 뚝 잘라서 FFT를 수행하면, 그 잘린 경계면의 영향으로 주파수 영역에서는 바늘처럼 뾰족한 단일 선 하나가 아니라 양옆으로 물결치며 퍼지는 모양(Sinc 함수)이 나타납니다. 수학적으로 표현하면 원래 신호에 사각 윈도우(Rectangular Window)를 곱해 잘라낸 것이 주파수 도메인에서는 컨볼루션(합성곱)이 되어, 윈도우의 주파수 응답인 Sinc 함수(Dirichlet Kernel) 형태로 퍼지게 되는 것입니다.

• 주파수가 정수 Bin에 딱 떨어질 때 ($k_0 = 8.0\text{MHz}$):

  Sinc 함수의 중심 피크가 정확히 $8$번 Bin에 오고, 다른 모든 정수 Bin($7, 6, 9, 10\dots$)들은 Sinc 함수의 영점(Zero-crossing: 값이 0이 되는 지점)에 정확히 걸치게 됩니다. 그래서 8번 Bin에만 에너지가 100% 집중되고 나머지 Bin은 완전히 0이 됩니다.

• 주파수가 Bin 사이에 걸칠 때 ($k_0 = 8.5\text{MHz}$):

  Sinc 곡선 자체가 옆으로 이동하면서 정수 Bin들이 영점을 벗어나 Sinc의 옆구리(Sidelobe)와 경사면을 샘플링합니다. 그 결과 피크 진폭은 약 0.64로 줄어들고(최대 -3.92 dB의 Scalloping Loss), 에너지가 인접한 여러 Bin으로 흩뿌려집니다. 이것이 바로 스펙트럼 누설(Spectral Leakage)입니다.

제미나이를 이용해서 아래 시뮬레이터를 만들었습니다. 슬라이더를 움직이거나 주파수를 입력하면, 해당 sine wave가 FFT를 거치면서 어떻게 인접한 bin에 영향을 미치는 지 직관적으로 이해할 수 있습니다.

이 원리를 알고 나면 반대의 질문, 즉 "작은 FFT 블록 사이즈로도 어떻게 1Hz 단위의 정밀한 주파수를 만들어낼 수 있는가?"에 대한 답이 명확해집니다. 굳이 5천만 개의 거대한 Bin을 만들지 않더라도, 인접한 몇 개의 Bin에 적절한 크기(진폭)와 위상 가중치를 배분하여 역변환(iFFT)을 수행하면 Bin 사이의 임의의 아날로그 주파수를 얼마든지 시간 도메인에 합성해낼 수 있습니다.

신호 정렬 문제

그런데 여기서 한가지 문제점이 생깁니다. 실제로 iFFT를 사용해 딱 떨어지지 않는 주파수를 생성해 보면, 프레임 사이에 불연속 구간이 생기는 것을 보게 됩니다.

FFT 프레임 경계면 불연속과 스펙트럼 누설

1. 정수 주기 신호 (A: Aligned, $k = 3.0$):

   프레임 1의 끝점($n=31$)과 다음 프레임 2의 시작점($n=32$)이 기울기와 값 모두 매끄럽게 이어집니다(연속). 따라서 주기적으로 무한히 이어붙여도 아무런 왜곡이 없고, 오직 하나의 주파수 성분만 검출됩니다.

2. 비정수 주기 신호 (B: Unaligned, $k = 3.5$):

   프레임 안에서 주기가 끝나지 못하고 골($-1.0$)에서 잘린 채 끝나는데, 다음 프레임은 다시 처음 시작값인 마루($+1.0$)부터 반복됩니다. 그 결과 프레임이 맞닿는 경계면마다 $\Delta y = 2.0$에 달하는 급격한 스텝 불연속이 강제로 발생합니다.

시간 도메인에서 이렇게 칼로 자른 듯한 단절은 푸리에 변환 시 임펄스/사각파와 같은 광대역 고주파 성분을 뿜어냅니다. 우리가 원하는 파형은 아니죠. 한 프레임 안에서는 주파수가 맞지만 전체적으로 보면 저 단절때문에 고주파 노이즈가 생기게 됩니다.

그러면 간단히 다음 프레임에서 파형을 생성할 때 이전 프레임에서 끝나는 지점에서부터 만들면 되겠죠? 이건 어떻게 구현할 수 있을까요?

이것도 실제로 해당 파형을 FFT 변환을 해 보면 알 수 있습니다.

우리가 원하는 것은 파형의 시작점($x[0]$)을 항상 $1.0$으로 고정하는 것이 아니라, 이전 프레임이 끝난 지점의 값에 맞춰 임의의 초기 위상($\phi$)을 부여하는 것입니다:

$$x[n] = \cos\left(2\pi \frac{k}{N} n + \phi\right)$$

이 파형을 FFT 변환하면 해당 주파수 Bin($k$)에는 단순한 진폭이 아닌 복소수(Complex Number: $\text{Re} + j\text{Im}$)가 들어갑니다. 오일러 공식에 의해 위상 $\phi$는 복소 평면 위에서 회전하는 벡터의 각도가 됩니다:

$$X[k] = A e^{j\phi} = A(\cos\phi + j\sin\phi)$$

• $\phi = 0^\circ$ (시작점 $x[0] = +1.0$): 실수부만 $+1.0$이고 허수부는 $0$입니다 (순수 코사인파).
• $\phi = -90^\circ$ (시작점 $x[0] = 0.0$, 상승): 실수부는 $0$이 되고 에너지가 허수부($+1.0j$)로 전부 이동합니다 (순수 사인파).
• $\phi = 180^\circ$ (시작점 $x[0] = -1.0$): 실수부가 $-1.0$, 허수부는 $0$이 됩니다 (반전 코사인파).
• 임의의 각도 $\phi$: 신호의 전체 크기(진폭 $|X| = \sqrt{\text{Re}^2 + \text{Im}^2} = 1.0$)는 절대 변하지 않고, 오직 실수부와 허수부의 비율(배분)만 달라집니다.

아래 시뮬레이터에서 위상(Phase) 슬라이더를 움직여보세요. 파형의 시작점이 바뀜에 따라 복소 평면에서 벡터가 어떻게 회전하고, FFT Bin의 실수부(파란색)와 허수부(주황색) 기둥이 어떻게 교차하는지 직접 확인할 수 있습니다.

결국 하드웨어나 소프트웨어로 연속적인 임의 주파수를 스트리밍 생성할 때의 해법이 여기에 있습니다.

"이전 프레임이 끝난 시점의 위상각 $\phi_{end}$를 누적해 두었다가, 다음 프레임의 FFT 주파수 Bin에 복소수 회전 계수 $e^{j\phi_{end}}$를 곱해서 iFFT를 수행하는 것"입니다. 복소수의 위상 회전만 적용하면 비정수 주파수라도 프레임 경계에서 뚝 끊기는 불연속점 없이 부드럽게 이어지는 연속 신호를 만들어낼 수 있습니다.

하드웨어 구현의 한계

이론상으로는 비정수 주파수는 복소수로 위상을 계속 추적해 주기만 하면 프레임의 경계면이 불연속되지 않고 연속적인 신호를 만들어낼 수 있습니다. 하지만, 실제 하드웨어나 소프트웨어로 구현할 때는 여러가지 한계가 발생합니다.

문제는 앞서 본 Sinc 형태로 퍼진 값을 각 Bin에 역으로 뿌려주는 과정(Sinc Splatting)입니다. 비정수 주파수의 스펙트럼 누설은 전체 FFT Bin에 영향을 미칩니다. 하지만 이걸 하드웨어로 구현하려면 모든 bin에 다 더해주어야 합니다. 이걸 주파수 하나만 생성하는 데에는 연산량이 얼마 되지 않습니다. 그러나 여러개의 주파수를 동시에 생성하는 경우에는 모든 Bin이 모든 주파수만큼 덧셈을 해야 합니다.

위의 그래프에서 보다시피, 비정수 주파수가 모든 Bin에 영향을 미치긴 하지만, 크게 영향을 받는 것은 인근 몇개의 Bin 뿐입니다. 그래서 하드웨어 연산량을 줄이기 위해 인근 4개 ~ 8개의 Bin만 더해주는 방식을 사용합니다. 요구하는 Signal-to-Noise (SNR) 에 따라 달라지겠죠.

이렇게 한정된 갯수의 Bin만 더하면 어떤 문제가 생길까요? 노이즈가 생기게 되고 프레임이 끝나는 부분에서 여전히 불연속점이 발생합니다.

이를 해결하기 위해 프레임을 겹쳐서 가중치를 주어 더하는 WOLA(Weighted Overlap-Add) 기법을 사용합니다. 겹치는 구간에서 이전 프레임 값은 점점 줄어들고(Fade-out), 다음 프레임 값은 점점 커지게(Fade-in) 만들어서, 전체 파형이 부드럽게 이어지도록 만듭니다. 여전히 근사치를 계산하기에 완벽하지는 않지만, 충분히 좋은 성능을 내 줍니다.

WOLA(Weighted Overlap-Add) 중첩 합산 기법

위 그림처럼 연속한 두 프레임을 50% 중첩(Overlap)시키고, 중첩 구간에서 이전 프레임에는 감쇄 윈도우($1 \to 0$)를, 다음 프레임에는 증가 윈도우($0 \to 1$)를 적용한 뒤 더합니다(Overlap-Add). 두 윈도우의 가중치 합이 항상 1.0으로 유지되는 조건(COLA: Constant Overlap-Add)을 만족하면, 잘려나간 Sinc 꼬리로 인한 프레임 경계 단절과 스텝 잡음을 없애고 연속적이고 깨끗한 파형을 합성할 수 있습니다.

마치며

대학교 시절 2D 이미지의 노이즈를 지우기 위해 처음 배웠던 FFT는, 시간이 흘러 하드웨어에서 실시간으로 연속 파형을 만들어내는 전혀 다른 문제로 다시 저를 찾아왔습니다.

처음에는 "원하는 주파수 해상도만큼 무조건 FFT 블록 사이즈를 키워야 하나?"라는 단순한 의문에서 출발했지만, 실제로 하드웨어 구현을 고민하며 파고들다 보니 **스펙트럼 누설(Sinc Splatting)**을 역이용해 Bin 사이의 정밀한 주파수를 표현하고, 프레임 경계의 불연속을 막기 위해 **복소 평면에서의 위상 누적($e^{j\phi}$)**을 추적하며, 하드웨어 리소스 한계로 잘려나간 Sinc 꼬리를 **WOLA(Weighted Overlap-Add)**로 부드럽게 봉합하는 과정까지 이어졌습니다.

수식으로만 보면 복잡하고 추상적인 푸리에 변환이지만, 결국 엔지니어링의 관점에서는 "한정된 하드웨어 자원(제한된 Bin 개수와 연산기) 안에서 어떻게 끊김 없고 깨끗한 연속 신호를 근사해 낼 것인가" 라는 트레이드오프였습니다. 예전에 잃어버렸던 블로그 글을 다시 쓰기 위해 시작했지만, 덕분에 교과서 속의 수식을 넘어 실제 실리콘 위에서 동작하는 FFT를 한층 더 깊이 이해할 수 있는 계기가 되었습니다.