디지털 신호처리

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

2. chapter 06. Fast Fourier Transform

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

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

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

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

7. 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에 대한 값은 자동적으로 알게 된다. ∴ 계산량은 현저히 감소함.

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

9. 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로 분해한 선도

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

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

12. 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]

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

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

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

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

17. 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로 분할함.

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

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

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

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

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

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

24. 6.6 FFT 응용  10[Hz]에서 스펙트럼의 최대치인 0[dB]이 되지 않고 1[Hz]부근에 상당한 피크치가 출현, 트랜드를 포함하지 않는 신호의 스펙트럼은 10[Hz]에서만 피크치가 존재함.

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

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

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

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

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

30. 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]부근의 대역에서 가장 큰 값을 갖는다.

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

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

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

34. 6.6 FFT 응용 자기상관함수의 계산은 N의 값이 클 때는 많은 시간이 소요 FFT를 이용하면 계산시간을 상당히 단축할 수 있다. FFT를 이용한 자기상관함수를 구하는 순서(그림 6.18) ⇒단, FFT를 사용하는 경우에는 N개의 데이터가 주기함수라야 함. 그림 6.18 자기상관함수 계산의 순서

35. 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문자를 역순으로 한 조어(造語). ⇒켑스트럼 해석 : 지진파, 음성해석 혹은 영상처리 등에 사용.