[경고] 아래 글을 읽지 않고 "고속 푸리에 변환"을 보면 바보로 느껴질 수 있습니다.
[고속 푸리에 변환의 쉬운 이해]
통신 이론과 신호 처리에 필수적인 고속 푸리에 변환(fast Fourier transform, FFT)은 이산 푸리에 변환(discrete Fourier transform, DFT)의 항들을 우함수(even function)와 기함수(odd function)로 모아서 급수를 빠르게 계산하기 위한 혁신적인 알고리즘이다[1], [2]. FFT는 1805년가우스 28세, 조선 순조 시절 무렵에 가우스Carl Friedrich Gauss(1777–1855)가 먼저 고안해 사용했다는 전설을 가지고 있다.[푸리에 교수가 푸리에 급수(Fourier series)를 공개한 시점이 1807년이다.] 현대적인 FFT는 뼈대 있는 수학자 룽에Carl Runge(1856–1927)의 1903년룽에 47세, 대한제국 시절 기여를 거쳐 미국 수학자 쿨리James Cooley(1926–2016)와 투키John Tukey(1915–2000)가 1965년쿨리 39세, 박정희 정부 시절에 재발견했다[3]. 쿨리와 투키가 FFT를 구현한 방법론은 쿨리–투키 알고리즘(Cooley–Tukey algorithm)으로 부른다. 이 알고리즘의 유도는 기존 DFT에 대한 관찰로부터 시작한다.
![](https://blogger.googleusercontent.com/img/b/R29vZ2xl/AVvXsEi2e3EwaX-mQDw7egpk3UItaPDFRGPRkWi9ekT6W8wjHyiN25PPmDDHXe2VLKm8vACNw6JvsTgK6jcUOx-XhP5IYW8O_M9rjv5uvv3Yr33DMBKRv6AT_USBJkzK3pBPENFhQifpS0jICxYu/s16000/dft5.png)
여기서 $M$은 표본 개수, $f_m$은 시계열(time sequence), $F_k$는 푸리에 계수(Fourier coefficient)이다. DFT의 모든 항을 항상 이등분할 수 있도록, FFT에서는 보통 $M$을 2의 거듭제곱인 $M$ = $2^N$으로 가정한다. 여기서 $N$은 1이상의 자연수이다. 먼저 식 (1)에서 복소 지수 함수를 $W_M$으로 치환해서 멱급수(power series) 형태로 다시 쓴다.
![](https://blogger.googleusercontent.com/img/b/R29vZ2xl/AVvXsEhORTZlkI7UDipTos1q0N82kgj0Ai5u1Syq760OZP4dBrhkdrQBvuS5KVQD89XaeiOAcGV90CW83tmMU2iMAAsMkrpi1ippcHwGIVwnDrH--fxJdtEAT7gHmlJmnCA-B66vXOuPwJIvZnV8MGTDmlNv1mF4jj3aaXog7Zn7yEGRradBTTIEZQawE2pMMw/s16000/ff1.png)
복소 지수 함수 $W_M$의 거듭제곱은 다음 성질을 가진다.
![](https://blogger.googleusercontent.com/img/b/R29vZ2xl/AVvXsEhRBDJXDUuEnlLkKUfS_sqGjG_fMgOMIsAXoRy5EAaRl5RY4uolYLqh4Rk3Zu5pBHj0qF8q9KIp6jnZL6kecq-WYRQf01p9fAqLK6Aya8wJb0FZP-5PyuCkn9AYBOu_s9Cok39CRdd8m86zp7esZLWRUCJkcpMhPlXhrptOlZaCugqV69i8Nm6uLzJyMg/s16000/ff2.png)
![](https://blogger.googleusercontent.com/img/b/R29vZ2xl/AVvXsEguBHHDkJ_3IxXTnuvbAayhQ4O11FxwFDx2QK5SqBV2LbdZgxP_GnV54SN9zuC-vESc4bpK493kpWlL0tDAn4pmTg540zU4OIivq3RBcugtWIj8BaMbPUYxhhyLyClfeGWwv-GCrrmp6-CZb12U5n48AfACZy8wRYX2-mMNyA-EpIy-V_jCLiQ5EncgBg/s16000/ff13.png)
그 다음에 식 (2)를 이등분해서 DFT 항을 $W_M^k$에 대한 우함수와 기함수인 부분으로 나눈다.
![](https://blogger.googleusercontent.com/img/b/R29vZ2xl/AVvXsEgCGKSkrVOrGdyn7BCALlX_lNB9upX_HoL4GSee3cf_nQozTH4Na_mVI79oscNg4_WhbvlTUWpRGmbgx058ViMqnjS_BpMhBAd2tE4TleySGr6Kysh6ZGNMRAiNljxomh7lvbTgk1oRKB2G2pLIaT-KTLa1RSHHPG1q_pHAhHGA0ODPAYk5b7muugcRKQ/s16000/ff14.png)
![](https://blogger.googleusercontent.com/img/b/R29vZ2xl/AVvXsEgqNEUtzIpvoDVS78_krHERv_6N-AGWuJ1xmwhYTK0SXkO7_vEG62nhsJLqrqcfTTxYgc1jVWgW0-wOH04KT7976QR5vmLoxa0OirQsAZF39001Hs0lItNcBqRVkI38tRbFcSFR9s4yKdsuOHPCYkgrQJLbioNS1BpYZiu7kY_V0Oqvt_UUIseaCeY92w/s16000/ff4.png)
여기서 $k$ = $0, 1, \cdots, M/2-1$, 원래 DFT 급수를 1번 분해한 $E_1(x)$와 $O_1(x)$는 각각 짝수 및 홀수 밑수를 가진 $f_k$만을 모은 다항식이다. 그러면 식 (3a)의 둘째식을 이용해서 $F_{k+M/2}$를 다시 쓸 수 있다.
![](https://blogger.googleusercontent.com/img/b/R29vZ2xl/AVvXsEjuZEB5WoKyaO--d_jQ5ow1zpjjf7VKzQZNSQIow7AlzdhGD-53T1MBAqXyeB_jFEMS8MoGqZoR81tIytRxILqujS4hBaAEJc6CjjqJ4tXCcSs4gYU8Y858VVErh-Xe0Dusxcw4-JZ6wK3-0vMgaZrF_mLBHzWv1P48k6-9RePeNKzAjzkJt09Rwg5VAQ/s16000/ff15.png)
여기서 $W_M^{k+M/2}$ = $-W_M^k$이다. 그러면 식 (6)을 다시 계산할 필요없이, 식 (4)에서 구한 $E_1(\cdot)$과 $O_1(\cdot)$의 결과값에서 부호만 ($+$)에서 ($-$)로 바꾸면 된다. 즉, DFT 급수를 우함수와 기함수로 분해함으로써 계산 속도를 2배만큼 증가시킬 수 있다. 이 방식을 염두에 두고 식 (4)를 새롭게 보면, $E_1(x)$와 $O_1(x)$는 항의 개수가 반[= $M/2$]만큼 줄어든 또 다른 DFT 급수이다. 그래서 재귀적인 방법으로 식 (5)를 다시 우함수와 기함수 관련 다항식으로 모아서 계속 연산을 가속화시킨다. 이와 같은 기법이 바로 쿨리–투키 알고리즘의 핵심이. 예를 들어, 재귀 관계에서 $E_1(\cdot)$과 $O_1(\cdot)$은 다음 다항식으로 재분해된다.
![](https://blogger.googleusercontent.com/img/b/R29vZ2xl/AVvXsEh6O9DrAeDnxTiixFLWvAzxbBprISYhfnSammM0KLRsRkokiZNOmg-OLqq8d76uCGCxILM1RilP_T3u6MobnggPaUp1i8HMvtiikrxxWS6Yc8ba0T7iVFzwlVXbLenodgOdsd-kzuqWRCMC-4M_qO-6k89HpbQiRWmwWre3CpaWFw2T4bMSaUZ_zbVILw/s16000/ff16.png)
![](https://blogger.googleusercontent.com/img/b/R29vZ2xl/AVvXsEguLBlvUU_5nLvN0LJ2iVC97eiukLh5yDeMR4QnMdVrKPuayeDV_8OzZXNoesMno3YXnyuEP5dsM2vo79T4yv1PsY26ggpRyzysRZM0k2ZqBY_7GRqA1-6ud7Vk1UX64FEZQDZyIcJnBdTYrRL6G66Z3P5ISB1loELUNET8YYTt0WKJ0XJR7giXND7XQg/s16000/ff17.png)
![](https://blogger.googleusercontent.com/img/b/R29vZ2xl/AVvXsEiqjj5iMgCn99xmiXO9sRUK-iw6-YLjhaU_ngP2O-VqNtzBc7bxnfBrMZdVdlwkI_qiidptfDFoK9aRSSsVspRx5jc5C0h84d3t7FFefywlsDPfRvjwXS-nVUMcAE2wz6DguHOug4A0hee_R_MIOxCB1ESX1GMZK5cT2lNXFcc7UAMW7C_pisiIHfi-6g/s16000/ff18.png)
![](https://blogger.googleusercontent.com/img/b/R29vZ2xl/AVvXsEiGIWxdXy8t7srjYzerVhRJbr_CYXG8PG50j0vkuqeMFNj83gOBkAz9dRpKmSfxS24nFE0Rbi93sgkhQ-KWnNt34NUDiu_t0jd6jFvL0wpIdY6U5-w1LfgiyhjfC9fggdGv_2XghYXZtgCPVCzif2f8ywM4c-f61xMvsMmEGOS5Md3hHQBEwIoeVOdFsg/s16000/ff19.png)
여기서 $k$ = $0, 1, \cdots, M/4-1$이다. 식 (8a)에서 $E_2 | E_1 (x)$와 $O_2 | E_1 (x)$는 각각 $E_1(x)$를 구성하는 $x$의 짝수차와 홀수차 항을 모아서 만든 다항식이다. 비슷하게 $E_2 | O_1 (x)$와 $O_2 | O_1 (x)$는 각각 $O_1(x)$의 짝수차와 홀수차 항으로 구성한다.
좀더 쉽게 쿨리–투키 알고리즘을 이해하기 위해, 예시적으로 $M$ = $2^1$인 경우부터 출발해 식 (4)를 다시 기술한다.
![](https://blogger.googleusercontent.com/img/b/R29vZ2xl/AVvXsEhnlxom7D3L4DwS9k2ykWoKYLVyju5U14UFa8BSDbDNod9GNcr3gQpIV7S6yUWItD5L5-NjmtNlxwHEz8bdgdbbm28ggh9d5DBnbBtDHICli8UvMjpH1W4VT0CScalsZd_V6m9jpkyUd7So-qdTQPiLZ30OJMbaQGTA98unym_GKg2zx3ZpHeMWZXAyqw/s16000/ff6.png)
![](https://blogger.googleusercontent.com/img/b/R29vZ2xl/AVvXsEhTmQZTdoia5hbE3PwObLLhW0HqYkJQYQTzBC-j44y-CvRqMKa6zXDhxSq8KeVahENb-cLoPNSLHnrY6AfPztJbRzKFfLoR9TKnJDaELzn9xPAV5zdKO_Lt1rojcecKWxbuaDHoLKAryyj58lOaG9-OipySp9yakmf0WrhLUaBlZxPmGT57AHeTydm1Qg/s16000/ff10.png)
[그림 1] 기본 $2 \times 2$ 나비 도표(사진 출처: wikipedia.org)
![](https://blogger.googleusercontent.com/img/b/R29vZ2xl/AVvXsEhoh5NRsmFKlIveNQeX95vrugSbudNAOAol7k9aV5aWpy6ndFx62yR3FzDM7vQJP--WYYa9BKWth4R6o6WBNXajfQ92COQeSNsG3CXc234mry4vDhWkQn3jgFELKDPeuhQ7pwBPx7H03_fhpz8Sq0SUTIzTqSz9F1wvhVOHigWdLCoNobdEZrgMSEQgdg/s16000/ff11.png)
[그림 2] 변형 $2 \times 2$ 나비 도표
[그림 1]은 나비를 닮아서 나비 도표(butterfly diagram)라고 명한다. 기본 나비 도표의 $f_1$에 $W$를 곱한 변형은 [그림 2]처럼 화살표를 하나 추가해서 그린다. 표본 크기를 더 증가시킨 $M$ = $4$인 경우도 식 (4)를 바탕으로 우함수와 기함수를 분리한다.
![](https://blogger.googleusercontent.com/img/b/R29vZ2xl/AVvXsEgLi2DsbZA_hh-2BXzrGvhkG5h3FeF6j7dYrXxX47FaDNuMvDbEQKI5GL1RLClaDjZMqiEXcjgF8jS7sQ96HEDeqCQzP88bHoGyfIRU7VmIBBkzP4QYWAvxuPMZDWV1Pof_j4J_Mis8Ot3ijiRb2i9zj32tdCNZw72NLWl3tRrUQBeO2BeA9iMynKywKA/s16000/ff9.png)
![](https://blogger.googleusercontent.com/img/b/R29vZ2xl/AVvXsEg1gjVBEftB84biX9inVYfcmVs1uDEtLY99VM-BlR40rjrBYeYqIASgIxoZZ2IbYa9uDMOeac8KE9QN7OZaFccjLzFrc08_y33Of1a4f3sa_op33A9Yp_Y9lzj_YJGG_Dy_V5a4rhaT2gUv0E9XDBH9jjYExsgK8zt8NayDR5TejwSjlL-cGHIGxWugMQ/s16000/ff12.png)
[그림 3] $4 \times 4$ 나비 도표
식 (10)을 나비 도표로 나타낸 결과인 [그림 3]은 약간 복잡하지만 명확한 규칙이 보인다. 시계열 $f_k$를 시작할 때, $k$를 이진수로 나타내어 비트를 역전시킨다. 이 과정은 FFT의 비트 역전(bit reversal)이라 이름 붙인다. 그러면 [그림 3]의 입력은 0, 2, 1, 3의 순서로 들어간다. 그 다음에 [그림 1]처럼 $M$ = $2$에 대한 FFT를 각각 수행하고, 이 결과를 모아서 각 화살표에 대해 $M$ = $4$에 대한 FFT 연산을 다시 적용한다. 식 (4)와 (7)은 비트 역전이 필요한 이유를 잘 보여주고 있다. 우리가 홀짝을 나눈 다음에 표본 개수 $M$과 지수 $k$를 1/2하고, 우함수와 기함수 관련 다항식을 각각 홀짝으로 다시 분류한다. 지수 $k$를 줄이는 단계는 2의 거듭제곱으로 나누면서 계속 반복된다. 최종적으로 항이 2개만 남으면 멈춘 후, 나비 도표인 [그림 1]을 이용해 계산한다. 시계열 $f_k$의 $k$가 줄어서 마지막에 선택된 두 항[그림 1에서 0과 1로 표기]은 $k$를 이진수로 표시할 때 최상위 비트(most significant bit, MSB)가 0이나 1인 경우이다. [그림 1]을 포함한 [그림 3]에서 나비 도표의 둘째 단계는 MSB의 그 다음 비트가 0이나 1을 가진 $k$를 처리한다. 따라서 FFT에서는 [그림 3]처럼 비트 역전을 통해서 입력 시계열 $f_k$를 두 개씩 모아야 한다.
![](https://blogger.googleusercontent.com/img/b/R29vZ2xl/AVvXsEjO4_67ocgY57qnS--4hv1jzQ0tkHgaZe-1SUgHDyhBMcrumasAHKhxJhshC4nqZljuWKfkZ0WdGJ-vF57GZrSWvbEdqk0Q2laln7doeupUgqOGWa3Dc1Whoh8SdqUVg1sBU1ELOY52IA71xopt9YQhRvhx5H4XWdEXGMj7JixD_lYBrD1SIXIkmkwkaQ/s16000/ff1.png)
[그림 4] $8 \times 8$ 나비 도표
[그림 4]는 [그림 3]을 2개 나열해서 $8 \times 8$ 나비 도표를 구현한 결과를 보여준다. 시계열 $f_k$의 $k$는 비트 역전을 통해서 결정한다. 입력이 8개라서 우함수와 기함수 부분인 다항식 $E_1(\cdot)$과 $O_1(\cdot)$은 각각 4개의 표본을 가지고 있다. 식 (4)를 사용해서 유도한 $M$ = $8$인 FFT 관계식은 다음과 같다.
![](https://blogger.googleusercontent.com/img/b/R29vZ2xl/AVvXsEgiGC4X_DmhLoISi5RzV_luX_qwkLm4tE28Ewn0XN95Wy6KokZtbCn34kSjl3RHTrGCx1PHxSafixwQ4NcRmthwXD1LLNEyKUcq9YFzaSJBy6H94VJxdt82pyvCwGS4f18LA6XPEYjbVCWgqtY5UqpH-LgZ9QNAFJhmaj1ARn0uqMmceU09UO11N3TpPw/s16000/ff2.png)
여기서 $k$ = $0, 1, 2, 3$, $E_1(\cdot)$과 $O_1(\cdot)$은 [그림 3]으로 계산한다.
표본 개수가 $M$ = $2^N$이 아니더라도 [1]을 참고해서 FFT를 수행할 수 있다. 하지만 $M$이 2의 거듭제곱이 아니면 홀짝 모으기가 번거로워진다. 만약 표본이 작아서 $M$ < $2^N$이면, $M$을 증가시켜 강제로 $M$ = $2^N$을 만들고 모르는 $f_k$는 $0$으로 채워서 FFT를 돌린다. 이런 부가적인 과정은 FFT의 영 채우기(zero padding)로 이름 붙인다. 영 채우기의 조건은 붙여넣은 $f_k$의 크기가 매우 작아서 0으로 근사할 수 있다는 가정이다. 예를 들어, $M_0$개로 영 채우기를 하면 표본수는 $M+M_0$으로 늘고 전체 표본화 시간은 $T$ = $(M+M_0) \Delta t$로 증가하다. 여기서 $\Delta t$는 표본화 주기(sampling period)이다. 이로 인해 주파수 영역의 각주파수는 $p_k$ = $k \Delta p$로 변한다. 여기서 $k$ = $0,1$, $\cdots$, $M+M_0-2, M+M_0-1$; $\Delta p$ = $2 \pi \mathbin{/} [(M+M_0) \Delta t]$가 되어 이전보다 줄어든다. 그래서 시간 영역의 표본수를 늘리면, 그만큼 주파수 영역의 표본 간격이 축소된다. 하지만 $p_k$의 최대값은 $2 \pi / \Delta t - \Delta p$이므로, 영 채우기로는 $p_k$의 범위를 거의 늘릴 수 없다.[$\Delta p$가 줄어든 정도만 최대값이 약간 증가한다.] 주파수 영역에도 영 채우기가 가능하다. 표본화 주기 $\Delta t$를 고정하고 시간 영역과 동일한 논리를 적용할 때, $M_0$를 아무리 늘려도 $p_k$의 범위는 거의 늘지 않고 주파수 간격 $\Delta p$만 좁아진다. 시간 영역에서는 전체 표본화 시간이 영 채우기한 만큼 늘어난다.
FFT의 기본 알고리즘은 간단해서 우리가 쉽게 프로그램할 수 있다. 그런데 많은 자료를 진짜 빠르게 처리하려면 전문적인 FFT 함수 집합 혹은 라이브러리(library)를 써야 한다. 우리가 쉽게 구할 수 있는 상용 FFT 라이브러리 중에 FFTW(FFT in the West)가 있다[4]. FFTW는 MATLAB에도 쓰일 정도로 잘 검증된 연산을 제공한다. 다만 FFTW는 비상업적인 용도로 쓸 때에만 무료로 이용할 수 있다. NVIDIA GPU를 사용하는 경우는 cuFFT[CUDA(Compute Unified Device Architecture) 기반 FFT]를 사용할 수 있다[5]. 여러 문헌에 나오는 결과를 종합하면, 표본수가 커질수록 cuFFT는 CPU를 쓰는 FFTW보다 훨씬 빠른 성능을 보인다.
[참고문헌]
[1] J. W. Cooley and J. W. Tukey, "An algorithm for the machine calculation of complex Fourier series," Math. Comput., vol. 19, no. 90, pp. 297–301, Apr. 1965.
[2] E. O. Brigham and R. E. Morrow, "The fast Fourier transform," IEEE Spectr., vol. 4, no. 12, pp. 63–70, Dec. 1967.
[3] M. T. Heideman, D. H. Johnson, and C. S. Burrus, "Gauss and the history of the fast Fourier transform," IEEE ASSP Mag., vol. 1, no. 4, pp. 14–21, Oct. 1984.
[4] M. Frigo and S. G. Johnson, "C Subroutine Library for Computing the DFT," FFTW 3.3.10, 2022. (방문일 2022-10-04)
[5] cuFFT Library User's Guide, cuFFT API Reference, NVIDIA Corp., California, USA, 2022. (방문일 2022-10-04)
댓글 없음 :
댓글 쓰기
욕설이나 스팸글은 삭제될 수 있습니다. [전파거북이]는 선플운동의 아름다운 인터넷을 지지합니다.