[Bio] #2 Matlab으로 DFT 구현하기

이미징 시스템을 처리하는 것은 결국 신호를 처리하는 것과 같다.

따라서 신호 처리에 기초가 되는

푸리에 변환(Fourier Transform)과 역변환(Inverse Fourier Transform)을 알아야 한다.

 

(신호 및 시스템 과목을 엉망으로 들었던 내겐 상당한 곤욕이었다.)

 

x(t)의 푸리에 변환

 

수학적 유도는 차치하고,

광학적인 관점에서 푸리에 변환이 어떻게 사용되는지 알아보자.

 

---

 

광학 분야 뿐 아니라 컴퓨터 비전(Computer Vision)에서는 특히,

푸리에 변환을 사용한다는 것이 

신호를 spatial domain에서 frequency domain으로 변환한다는 것을 의미한다.

 

 

 

Discrete한 domain에서 신호를 푸리에 변환할 때는 다음과 같은 공식을 사용할 수 있다.

 

DFT(Discrete Fourier Transform)

 

이를 이산 푸리에 변환(Discrete Fourier Transform, DFT)이라고 하는데,

이 이산 푸리에 변환과 그 역변환을 빠르게 수행하게 해주는 알고리즘을

'고속 푸리에 변환(Fast Fourier Transform, FFT)' 이라고 한다.

 

 

특정 알고리즘을 칭하는 말은 아니고,

DFT 연산을 빠르게 진행할 수 있는 알고리즘들을 지칭하는 말이라고 생각하면 된다.

대표적인 FFT에는 쿨리-튜키 알고리즘이 있다.

 

 

DFT 연산에는 O(n^2)이라는 시간이 걸리는데,

FFT를 사용하면 일반적으로 수행시간이 O(nlog(n))이 된다.

 

 

푸리에 변환은 영상 처리와 필터링을 할 때 필수적이다.

spatial domain의 신호를 frequency domain으로 변환하면

신호의 추가적인 특성들을 파악할 수 있기 때문이다.

 

---

Matlab Code

 

DFT(Discrete Fourier Transform)을 수행하는 matlab 코드를 작성해보자.

 

function G = DFT(g)
    M = length(g);
    W = exp(-1j*2*pi/M*(0:M-1)'*(0:M-1));
    
    G = W*g;
end

(티스토리 코드 블럭에 Matlab은 지원이 안돼서 따로 Javascript를 이용해 커스터마이징 해줘야 하는 듯하다.)

(그래서 엉뚱하게 하이라이팅 됐다.. 조만간 고쳐볼게요..)

 

 

 

코드는 무척 간단하다.

 

M*1 크기의 배열 g를 입력받으면(이때 type은 실수든 복소수든 상관없다.),

여기에 exp ( -1i * 2 * pi / M * (0:M-1)' * (0:M-1))을 곱해준다.

여기서 (0:M-1)' * (0:M-1)은 0부터 M-1까지 M개의 정수로 행, 열이 채워진 M*M 크기의 정수 행렬이다.

이 M*M 행렬의 각 요소에 -1i * 2 * pi / M 을 곱한 후 exp를 씌우고, 이를 다시 g와 곱해 Sigma를 구현했다.

 

그러면 g를 DFT한 결과인 G도 M*1 크기의 배열이 된다.

 

 

 

 

이 코드를 테스트 하기 위해 다음과 같은 과정을 거쳤다.

 

10*1 double matrix

 

 

시간을 잰 이유가 있다.

 

 

Matlab 내장 fft를 사용해 결과 검토

 

fft를 사용했을 때와 같은 결과를 나오는 걸 확인했다.

 

또한 내장 fft를 사용하니 경과 시간이 더 짧아진 것을 확인할 수 있다.

물론 내장 함수를 사용하면 최적화 등의 문제로

직접 구현한 함수보다 경과 시간이 느려지는 경우도 있는데,

내장 fft는 이름대로 내가 구현한 DFT보다 fast했다.

 

 

 

앞서 언급했듯이 FFT는 다양한 분야에서 사용된다.

앞으로 신호처리에 대해 정리할 일이 생기면 더 구체적으로 다뤄보도록 하겠다.