태그 보관물: Mathematics

Fast Fourier Transform

공학에서 FFT만큼 많이 나오는 용어도 흔치 않으리라 보입니다. 왜냐하면 공학에서 현상을 해석할 때에는 시간에 대한 해석 뿐만이 아니라 주파수에 대한 해석도 병행하기 때문입니다.(근래에는 두가지 영역을 동시에 해석하려는 시도도 나오고 있습니다.) 그런데 대다수의 사람들이 FFT를 많이들 어려워합니다. 본 강의는 FFT의 특성을 Matlab™ 을 통해 그래픽적으로 분석함으로서 정확한 FFT의 이해를 목적으로 합니다.먼저 강의를 진행하기 전에 일러둘 것이 있습니다. 강의에서 다룰 것은 실제적으로는 DFT(Discrete Fourier Transform)라는 것입니다.DFT는 FFT와 여러모로 유사하지만 다른 점이 있다면 FFT는 속도개선 차원의 알고리즘이 더 들어가 있다고 보시면 됩니다. FFT의 속도 향상에 관심이 있으신 분들은 다음 책과 사이트에서 도움을 받으실 수 있을 것입니다.

Richard G. Lyons 저 “Understanding Digital Signal Processing”
http://ocw.mit.edu/OcwWeb/Mathematics/18-06Linear-AlgebraFall2002/VideoLectures/index.htm  사이트의 26번째 비디오 강의

그 외에 DSP 전반의 시뮬레이션을 Matlab™ 을 통해 확인해보고 싶으신 분은 다음 두 책을 참조하시길 권합니다.

Proakis 저 “Digital Signal Processing Using Matlab”
Samuel D. Sterns 저 “Digital Signal Processing with Examples in Matlab”

먼저 FFT의 간단한 정의부터 확인해 보겠습니다.

Continuous Fourier Transform은 다음과 같이 정의됩니다.

그렇지만 디지털 개념으로 바뀌게 되면 다음과 같은 정의를 할 수가 있게 됩니다.


여기서 각각의 용어들은 다음과 같이 정의됩니다.

예를 들어보도록 하겠습니다.다음과 같은 신호가 있다가 합시다.만약에 N=8, fs=8000 samples/s 라고 한다면 Analysis Frequency는 m*fs/N = 1000*m 이 됩니다.결국 index에 1000 을 곱한 값이 되는 겁니다.

그럼 이 상태에서 FFT를 수행하도록 하겠습니다.

X(1)과 X(2)는 각각 1kHz와 2kHz이므로 값이 나와야 합니다. 그런데 X(6)과 X(7)에서도 성분이뜨는 것을 확인할 수가 있습니다. 뭔가 이상하다가 생각할 수도 있지만 FFT가 N/2 지점을 기준으로 대칭이라는 것만 아신다면 자연스러운 결과라는 것을 알 수가 있습니다.

그런데 위의 그림에서도 보여지듯이 원래는 1과 0.5의 크기가 떠야 할것이 4와 2의 크기로 생겼음을 알 수가 있습니다. 그것은 FFT가 원래 크기에 N/2를 곱한 크기를 나타내기 때문입니다.

자, 이쯤에서 Matlab™을 이용해 FFT를 직접 구현해 보도록 하겠습니다.Matlab Coding N=8;fs=8000;ts=1/fs;for n=1:1:N x(n)=sin(2*pi*1000*(n-1)*ts)+0.5*sin(2*pi*2000*(n-1)*ts+3*pi/4);endfor m=1:1:N temp=0; X=0; for n=1:1:N X=x(n)*(cos(2*pi*(n-1)*(m-1)/N)-j*sin(2*pi*(n-1)*(m-1)/N)); X=X+temp; temp=X; endXX(m)=X;endm=0:1:N-1;subplot(2,2,1)stem(m,abs(XX));subplot(2,2,2)stem(m,phase(XX)*180/pi);

Matlab Internal function(fft)

N=8;fs=8000;ts=1/fs;for n=1:1:N x(n)=sin(2*pi*1000*(n-1)*ts)+0.5*sin(2*pi*2000*(n-1)*ts+3*pi/4);endXX=fft(x);    %% or, XX=fft(x,N);m=0:1:N-1;subplot(2,2,3)stem(m,abs(XX));subplot(2,2,4)stem(m,phase(XX)*180/pi); 

위에서도 보여지듯이 Matlab™ 에서 제공하는 fft라는 함수를 이용하면 코드가 훨씬 간단해지는 것을 보실 수가 있습니다.(참고로 프로그램이 제대로 실행 안되시는 분들은 수정 파일을 첨부했으니 참조하셔서 실행해보시기 바랍니다.) 이제 본격적으로 FFT의 특성을 시뮬레이션을 통해서 확인해 보도록 하겠습니다. 신호와 조건은 다음 그림과 같습니다.


먼저 Sampling Frequency 의 중요성을 알아보도록 하겠습니다.
우선  Sampling Frequency 가 Nyquist Frequency 보다 작은 99일 때를 보도록 하겠습니다.처음 그림은

원래 신호이며 두번째 그림은 FFT한 그림, 세번째 그림은 FFT한 신호를 다시 IFFT를 통하여 되돌린그림입니다.

엉뚱하게도 50Hz 에 떠야 할 성분이 49Hz 에서 뜬 것을 확인할 수가 있습니다.

다음은 Sampling Frequency 가 Nyquist Frequency 와 같은 100일 때를 보도록 하겠습니다.

이번엔 50Hz의 성분 크기가 두배가 되어 버렸습니다. 그 원인은 다음 그림과 같이 전체적인 주파수 특성이 Nyquist Frequency/2 지점을 기준으로 대칭이기 때문입니다.


다음은 Sampling Frequency 가 Nyquist Frequency 보다 큰 101일 때를 보도록 하겠습니다.

올바른 값이 나온 것을 확인할 수가 있습니다.Sampling Frequency 가 102일 때도 다음 그림과 같이 올바른 값을 확인할 수가 있습니다.

이제는 주기 일치의 중요성을 알아보도록 하겠습니다.신호를 좀더 분명하게 보기 위하여 Sampling Frequency를 1000으로 놓고 시뮬레이션을 하도록 하겠습니다.이번에는 주기를 1주기,2주기,10주기까지 주어가며 확인해 보도록 하겠습니다.

1 Cycle

2 Cycles

10 Cycles

윗 그림에서  주기가 늘어날 수록 관심 주파수가 아닌 곳의 포인트 수가 점점 늘어나는 것을 확인해 볼 수가 있습니다.여기까지는 별다른 문제가 없어 보입니다. 그러나 주기가 4.4인 경우를 살펴보도록 하겠습니다.

그림에서도 보여지듯이 정확한 값을 얻는데 실패하였고 다른 성분들도 튀어나오는 것을 확인할 수가 있습니다. 참으로 난감한 일이 아닐 수 없습니다. 이렇게 주기를 일치 시켜줄 수 없을 때에  문제들을 해결하기 위한 몇가지 방법이 있습니다. 첫번째는 Windowing이라는 기법이며 두번째는 Zero padding 기법입니다. 먼저 Windowing 기법을 살펴보도록 하겠습니다. Windowing 은 다음 그림과 같은 신호 변형을 원신호에 가하게 됩니다.

Windowing 이 이번 예제에 쓰인 경우와 다른 신호에 쓰인 경우를 비교해보도록 하겠습니다.

CDMA2000

그림에서 볼 수 있듯이 기본적으로 Windowing은 CDMA2000 신호와 같이 대역폭을 지니는 신호에 적합합니다. 그러므로 본 예제에는 적합하지 않은 방식입니다.그러면 이제 Zero padding 기법을 보도록 하겠습니다.


기본적으로 위와 같은 방식이 Zero padding 입니다.Zero padding 은 Matlab™ 에서 fft(x,N);이라는 명령어로 구현할 수가 있습니다. 여기서 N에 Signal length보다 큰 값을 주게 되면 그게 바로 Zero padding 의 효과를 주게 됩니다.Zero padding 의 길이를 늘려줌에 따라 결과가 아래와 같이 나오게 됩니다.

N = 2*Signal length

N= 3*Signal length

N = 10*Signal length

주파수 성분은 잘 뽑아 내었지만 Zero padding 의 길이가 늘어남에 따라 필요없는 대역이 생겨났음을 확인할 수가 있습니다. 이럴 때는 다음 그림과 같이 주기를 늘려주면 올바른  결과를 얻게 됩니다.


결론적으로 올바른 fft를 위해서는 다음과 같은 규칙을 지켜줘야 한다는 것을 알 수가 있습니다.

1.Sampling Frequency는 Maximum Frequency보다 최소한 2배이상을 잡아주어야 합니다.
2.주기를 되도록이면 맞추어 주셔야 합니다.
3.주기를 많이 주실수록 좋습니다.
4.주기를 못 맞추었을때에는 Zero padding 을 이용하시면 됩니다.

[출처] Matlab을 이용한 FFT의 올바른 이해|작성자 일편단심