디지털 신호처리 Jhmoon93@gmail.com

Fast Fourier Transform chapter 06. Fast Fourier Transform

주요내용 개요 시간솎음 알고리즘 주파수솎음 알고리즘 IDFT 알고리즘 FFT 응용

6.1 개 요 고속푸리에변환(Fast Fourier transform : FFT) 6.1 개 요 고속푸리에변환(Fast Fourier transform : FFT)  이산푸리에 변환(DFT)을 고속으로 산출하기 위한 하나의 알고리즘.  시간솎음(시간영역 분해), 주파솎음(주파수영역분해)알고리즘 이산푸리에변환을 집적 계산할 경우, X(k)복소계산 → N·N=N2회의 곱셈과 N(N-1)회의 덧셈이 필요. N이 커질수록 계산량은 비약적으로 증가. 고속 푸리에 변환을 이용해 계산량을 산출하면 → (N log2N)/2회의 곱셈과 2log2N회의 덧셈이 필요, 계산시간은 곱셈의 회수에 좌우되므로 N의 값이 클수록 고속푸리에 변환의 필요성을 알수 있다.

6.1 개 요  데이터 수 N = 2r(r는 양의 정수)이라 할 때 DFT 직접계산과 FFT의 곱셈회수 비교

6.2 시간솎음 알고리즘  N점 DFT X(k)를 구할때, 수열 x(n)을 n이 짝수인 것과 여기서, N = 2r  식 (6.1)의 우변을 n이 짝수인 항과 홀수인 항으로 분류하면  식 (6.3)에서 우변의 제1,2항에 속하는 x(n)개수는 동일하며 n이 짝수일 때는 n = 2r, 홀수일 때는 n = 2r + 1로 표시.

6.2 시간솎음 알고리즘  X(k)는 N점 DFT이지만, G(k)와 H(k)는 모두 N/2점 DFT이다. 그런데 WN/2rk는 k에 관하여 주기가 N/2인 주기함수, 그 이유는  G(k)와 H(k)를 구할 때 k=0, 1,…, (N/2) - 1에 대한 값만 산출하면, k = N/2, (N/2)+1, … , N-1에 대한 값은 자동적으로 알게 된다. ∴ 계산량은 현저히 감소함.

6.2 시간솎음 알고리즘  식 (6.5)를 이용해 N=8일 때 X(k)를 구하는 과정의 흐름선도

6.2 시간솎음 알고리즘  g(r) = x(2r)이기 때문에 g(2 l) = x (4 l), g (2 l + 1) = x (4 l + 2)로 된다. 그림 6.2 N/2점 FFT로 분해한 선도

6.2 시간솎음 알고리즘 그림 6.3 그림6.2를 그림 6.1에 대입한 선도  WN/2 =WN2의 관계가 성립함.

6.2 시간솎음 알고리즘 그림 6.4 2점 FFT  N = 8일 경우의 계산량.  좌단에서 우단에 도착하는데 3단계 계산과정을 거침.  각 과정마다 WNk를 곱하는 회수는 8로 됨. ∴ N의 일반적인 값에 대해 좌단에서 우단에 도착하는 데는 log2 N단계를 거쳐야 하며 각 단계마다 WNk 를 곱하는회수는 N 그림 6.5 8점 FFT를 산출하는 과정

6.2 시간솎음 알고리즘 1. 역비트순 (bit-reversed order) 역비트순의 과정  DFT를 산출하는 과정에서 입력수열의 나열은 순서대로 되어 있지 않기 때문에 입력순서를 얻어내기 위한 방법. 역비트순의 과정 step 1 : x(n)의 n을 2진수로 표시한다. step 2 : bit의 순서를 역으로 한다. step 3 : 다시 10진수로 표시한다. 예) N = 8인 경우 8 = 2k 이므로 k = 3 [bit]

6.2 시간솎음 알고리즘 2. 나비선도 ∴  N이 주어졌을 때 DFT를 시간솎음 FFT로 산출할 경우, 각 계산단계에 따른 p , q , r 값의 결정 ⇒ 식 (6.6), (6.7)등 그림 6.7 기본적인 나비흐름 선도 그림 6.8 간소화된 나비흐름 선도  나비계산 (butterfly computation) ⇒ 여기서 ∴

6.2 시간솎음 알고리즘  일반적으로 각 단계마다 N/2개의 나비가 존재. 그림 6.9 나비계산을 이용한 FFT 흐름 선도  일반적으로 각 단계마다 N/2개의 나비가 존재.  하나의 나비에는 1회의 곱셈이 필요하므로 그림 6.9에 의한 DFT 산출에 필요한 곱셈의 총 회수는 (N/2) log2 N  DFT 산출에 필요한 총 곱셈회수 Nlog2N과 비교하면 곱셈회수가 절반으로 감소함.

6.3 주파수솎음 알고리즘 주파수솎음(decimation-in-frequency) FFT 알고리즘  출력 DFT X(k)를 순차적으로 산출할 수 있음.  DFT X(k)를 순차적으로 작은 부수열로 분할하는 것에 기초.  주파수라는 용어는 X(k)의 k가 주파수를 나타내기 때문임.  N = 2r 이라 하고, 입력수열 x(n)을 순서대로 위 절반과 아래 절반으로 나누면 (식 6.12) ⇒ 우변 제2항을 변수치환하면

6.3 주파수솎음 알고리즘 ⇒ 여기서, WN(N/2)k = (-1)k 이므로 식 (6.13)은  K = 0, 1, … , N-1이므로 식 (6.14)에서 X(k)를 k가 짝수인 것과 홀수인 것으로 나누어 쓰면 ⇒ WN2rn = WN/2rn의 관계에 의해

6.3 주파수솎음 알고리즘  주파수솎음에 의하여 8점 DFT를 2개의 4점 DFT로 분할할 경우, 식 (6.16) 및 식 (6.17)에 따라 DFT X(k)를 구하기 위한 DFT 흐름선도 그림 6.10 2개의 4점 FFT 흐름 선도  N점 DFT를 두 개의 N/2점 DFT로 분할함.

6.3 주파수솎음 알고리즘  N/2점 DFT를 2개의 N/4점 DFT로 분할할 수가 있다. ⇒ 여기에서 의 관계가 성립함. 그림 6.11 N/4점 FFT로 분할한 FFT 흐름 선도 ⇒ 여기에서 의 관계가 성립함.

6.3 주파수솎음 알고리즘  N = 2r의 경우 주파수솎음 알고리즘에 필요한 곱셈 ⇒ (N/2)log2N,  N/4점 DFT를 N/8점 DFT로 분할 그림 6.12 N=8의 주파수솎음 FFT  N = 2r의 경우 주파수솎음 알고리즘에 필요한 곱셈 ⇒ (N/2)log2N, 가산회수 ⇒ Nlog2N ∴ 총 계산량은 시간솎음의 경우와 동일함. 시간솎음 FFT 선도와 주파수솎음 FFT 선도는 전치(transpose)의 관계임.

6.4 IDFT 알고리즘  IFFT(Inverse Fast Fourier Transform) 알고리즘 ⇒ IDFT에 대한 x(n)의 산출 알고리즘.  DFT와 IDFT(Inverse DFT)의 식 그림 6.13 주파수솎음 IFFT 알고리즘

6.5 기타 FFT 알고리즘  FFT보다 더욱 고속인 WFT(Winograd Fourier Transform) 알고리즘. PFF(Prime Factor FFT) 알고리즘. ⇒ 데이터수가 서로 소(素)인 수의 곱으로 되어 있을 때 적용하는 알고리즘. ⇒ 승산회수가 종래 2를 기수로 하는 FFT보다 약 20∼30[%]로 줄어듬.  LSI 기술의 발달에 의해 승산의 연산이 가감산의 연산과 거의 같은 속도로 실행가능한 디지털 시그날 프로세서(DSP)가 개발.  DSP의 특징을 고려한 FFT 알고리즘 출시.  삼각함수대신 구형파 함수를 이용한 Walsh 변환.

6.6 FFT 응용 1. 스펙트럼분석 1) 저역통과 필터를 통과시킨다. 2) A/D 변환  어떠한 신호에 포함되어 있는 주파수성분의 분포를 구하는 것. 1) 저역통과 필터를 통과시킨다.  신호의 최고주파수를 맞추어 놓기 위하여 먼저 저역통과 필터를 통과 ⇒ 신호처리에서 대상이 되는 신호의 최고주파수를 사전에 정확히 아는 것은 상당히 어렵기 때문임. ∴ 저역통과 필터를 이용하여 필요한 정보가 없어지지 않을 정도로 대역제한을 하여 최고주파수를 강제적으로 결정함. 2) A/D 변환  컴퓨터에 신호를 입력하기 위하여 아날로그 신호를 디지털 신호로 변환.

6.6 FFT 응용 3) 평균값(신호의 직류성분)을 0으로 데이터 변환한다. 4) 트랜드 제거  디지털 신호의 데이터 수가 N일 경우 평균치를 구함.  각 표본치로부터 이 값을 뺀 값을 새로운 데이터로 설정함. 4) 트랜드 제거  트랜드(trend) : 해석할 데이터 길이보다 긴 주기를 가진 성분.  트랜드를 포함하게 되면 신호는 시간과 함께 변화함.

6.6 FFT 응용  10[Hz]에서 스펙트럼의 최대치인 0[dB]이 되지 않고 1[Hz] 부근에 상당한 피크치가 존재함.

6.6 FFT 응용 ∴ FFT에 의한 스펙트럼분석을 하기 전에 이것을 제거해야함.  Hanning 창함수를 이용해도 정확한 스펙트럼은 얻어지지 않는다. ⇒ 트랜드에 의해 스펙트럼에 커다란 오차가 발생하기 때문 ∴ FFT에 의한 스펙트럼분석을 하기 전에 이것을 제거해야함.

6.6 FFT 응용  일반적으로 트랜드를 제거하기 위해서 최소자승법을 많이 이용함.  표본화 간격 t[sec]로 A/D 변환한 N개의 데이터  데이터에 포함된 트랜드 trn이 x의 R차 다항식으로 표현.  최소자승법 ⇒ x(n)과 trn의 차(差)의 자승이 최소가 되도록 계수 ar를 구하는 것.  계수 ar 를 구하기 위해서

6.6 FFT 응용 선형 트랜드 성분을 제거하는 경우의 예 ⇒ 식 (6.21)에서 차수 R = 1의 경우 그림 6.15 최소자승법에 의한 트랜드 제거 ⇒ 식 (6.21)에서 차수 R = 1의 경우 ⇒ By 식 (6.23)으로부터 계수 a0 , a1 을 구하면

6.6 FFT 응용 그림 6.15(b) ⇒ 최소자승법에서 선형 트랜드 성분을 추정한 것. 그림 6.15 최소자승법에 의한 트랜드 제거 그림 6.15(b) ⇒ 최소자승법에서 선형 트랜드 성분을 추정한 것. ⇒ 계수 a0는 절편을, a1은 직선의 기울기를 나타냄. 그림 6.15(c) ⇒ 원 신호로부터 선형 트랜드를 제거한 것.

6.6 FFT 응용 5) 창함수 적용 6) FFT 처리  추출한 신호에 불연속점이발생할 가능성이 있을 경우에는 창함수를 이용하여 처리한다. 6) FFT 처리  전(前) 처리가 끝난 신호에 대하여 FFT에 의한 주파수분석이 가능.  샘플수 N을 2r으로 만들기 위해 나머지 부분에 0을 첨가시켜 보다 정확한 스펙트럼을 얻음. 예) N = 100인 경우 0을 28개 채워 샘플수를 128개로 함.

6.6 FFT 응용 그림 6.16 뇌파 스펙트럼 ① 성인이 안정하고 있을 때 뇌파의 스펙트럼.(그림 6.16) ② 0.1~30[Hz]의 대역통과 필터를 통과후, T=5[ms]의 표본화 간격으로 A/D변환 ③ FFT 처리의 샘플수는 256개(1.28[sec])로 함. ④ 뇌파는 불규칙 신호이므로 1.28[sec]의 파워스펙트럼을 구하여 50회 반복, 파워스펙트럼의 가산평균을 구함. ⑤ FFT의 해석시간은 1.28[sec]이므로 파워스펙트럼은 약 0.78[Hz] 간격마다 선스펙트럼이 나온다. ⑥ 안정 상태의 뇌파는10[Hz] 부근의 대역에서 가장 큰 값을 갖는다.

6.6 FFT 응용 2. 자기상관함수 계산 상호상관함수(cross correlation) 상관함수 자기상관함수(auto correlation)  상호상관함수  두 신호의 유사성, 시간차를 나타내며 잡음이 포함되어 있는 신호의 검출 및 복원 등에 응용. 예) 레이다(radar) 신호의 탐지, 패턴매칭(pattern matching), 지연측정 등에 사용  자기상관함수  상호상관함수의 특별한 형태로서 어떠한 신호가 잡음과 함께 있을 때 그 신호의 주기검출에 주로 이용

6.6 FFT 응용  자기상관함수는 시간적으로 연속인 신호에 대하여 정의됨. 그림 6.17 자기상관함수 그림 6.17 자기상관함수  자기상관함수는 시간적으로 연속인 신호에 대하여 정의됨. 여기서  ⇒ x(t)와 x(t+ )의 시간차

6.6 FFT 응용  표본화한 이산신호 x(n)(n = 0, 1, …, N-1)의 자기상관함수 여기서, m ⇒ x(n)이 x(n+m)과 얼마나 어긋나고 있나를 샘플점의 수로 표현한 것.  m에 표본화 간격 t를 곱하면 시간차가 된다.  자기상관함수의 스펙트럼 P XX(h)

6.6 FFT 응용  자기상관함수의 계산은 N의 값이 클 때는 많은 시간이 소요 주기함수라야 함. 그림 6.18 자기상관함수 계산의 순서

6.6 FFT 응용 3. 켑스트럼 계산  신호 x(n)의 켑스트(cepstrum) cx(m)  켑스트럼 X(k) : 길이 N인 이산신호 x(n)의 이산푸리에변환 m : 케프런시(quefrency) log|X(k)| : 대수진폭(log magnitude) 그림 6.19 켑스트럼 계산의 순서  켑스트럼 ⇒ 핏치의 영향을 없애고 평활한 스펙트럼 포락을 추출하는 방법. ⇒ 스펙트럼(spectrum)의 앞의 4문자를 역순으로 한 조어(造語). ⇒ 켑스트럼 해석 : 지진파, 음성해석 혹은 영상처리 등에 사용.