Main Content

푸리에 변환

푸리에 변환 정의

푸리에 변환은 영상을 다양한 크기, 주파수, 위상을 갖는 복소수 지수의 합으로 표현합니다. 푸리에 변환은 향상, 분석, 복원, 압축 등 다양한 영상 처리 분야에서 중요한 역할을 수행합니다.

f(m,n)이 2개의 이산 공간 변수 mn의 함수라면 f(m,n)2차원 푸리에 변환은 다음과 같은 관계로 정의됩니다.

F(ω1,ω2)=m=n=f(m,n)ejω1mejω2n.

변수 ω1ω2는 주파수 변수로, 단위는 샘플당 라디언입니다. F(ω12)f(m,n)주파수 영역 표현이라고 불리기도 합니다. F(ω12)ω1ω2가 모두 주기 2π를 갖는 복소수 값 함수입니다. 이러한 주기성 때문에 보통 범위 πω1,ω2π만 표시됩니다. F(0,0)f(m,n)의 모든 값의 합입니다. 이러한 이유로, F(0,0)은 푸리에 변환의 상수 성분 또는 DC 성분이라고 합니다. (DC는 직류를 나타내며, 전기공학에서 정현파 변전압을 갖는 전원이 아닌 정전압을 갖는 전원을 가리키는 용어입니다.)

변환의 역은 변환된 영상에 대해 수행했을 때 원본 영상이 생성되는 연산입니다. 2차원 푸리에 역변환은 다음과 같이 표현됩니다.

f(m,n)=14π2ω1=ππω2=ππF(ω1,ω2)ejω1mejω2ndω1dω2.

개략적으로 말해서 이 방정식은 f(m,n)을 서로 다른 주파수를 갖는 무한개의 복소수 지수(정현파)의 합으로 표현할 수 있음을 의미합니다. 주파수 12)에서의 비중의 크기와 위상은 F(ω12)로 표현됩니다.

푸리에 변환 시각화하기

이에 대한 예로, 사각형 영역 안에서는 값이 1이고 영역 밖에서는 0인 함수 f(m,n)이 있다고 가정하겠습니다. 도식을 단순화하기 위해, 변수 mn은 이산 변수이지만 f(m,n)은 연속 함수로 표시했습니다.

사각 함수

Plot of rectangular function f(m,n) in the spatial domain

다음 그림은 위 그림에 표시된 사각 함수의 아래와 같은 푸리에 변환 크기를 메시 플롯으로 보여줍니다.

|F(ω1,ω2)|,

크기를 메시 플롯으로 표시하는 것은 푸리에 변환을 시각화하는 일반적인 방법입니다.

사각 함수의 크기 영상

Mesh plot of the magnitude of the Fourier transform of the rectangular function f(m,n), plotted as a function of the horizontal frequency, ω1, and vertical frequency, ω2

플롯 중앙의 피크는 F(0,0)으로, 이것은 f(m,n)의 모든 값의 합입니다. 플롯을 통해 F(ω12)가 세로 고주파에서보다 가로 고주파에서 더 많은 에너지를 갖고 있다는 사실도 알 수 있습니다. 이는 f(m,n)의 가로 단면은 폭이 좁은 펄스이고 세로 단면은 폭이 넓은 펄스라는 사실을 반영합니다. 폭이 좁은 펄스는 폭이 넓은 펄스보다 더 많은 고주파 성분을 갖습니다.

푸리에 변환을 시각화하는 또 다른 일반적인 방법은 다음 그림과 같이 아래 식을 영상으로 표시하는 것입니다.

log|F(ω1,ω2)|

사각 함수의 푸리에 변환의 로그

2-D pot of the log of the Fourier transform of the rectangular function f(m,n), plotted as a function of the horizontal frequency, ω1, and vertical frequency, ω2

로그를 사용하면 F(ω12)가 0과 매우 가까운 영역에서 푸리에 변환의 세부 정보를 살펴보는 데 도움이 됩니다.

다음은 그 밖의 간단한 형태의 푸리에 변환 예입니다.

간단한 형태의 푸리에 변환

2-D plots of the Fourier transforms of a rectangle tilted forty-five degrees, a circle centered at (0, 0), and an X shape centered at (0, 0)

이산 푸리에 변환

일반적으로 컴퓨터에서 푸리에 변환을 수행할 때는 이산 푸리에 변환(DFT)이라는 형태의 변환을 사용하게 됩니다. 이산 변환은 입력값과 출력값이 이산 샘플이어서 컴퓨터가 조작하기 쉬운 변환입니다. 다음은 이러한 형태의 변환을 사용하는 두 가지 주된 이유입니다.

  • DFT의 입력값과 출력값은 모두 이산 값이므로 컴퓨터가 조작하기 쉽습니다.

  • DFT를 계산할 때 고속 푸리에 변환(FFT)이라는 빠른 알고리즘을 사용할 수 있습니다.

DFT는 일반적으로 유한 영역 0mM10nN1에서만 0이 아닌 이산 함수 f(m,n)에 대해 정의됩니다. 2차원 M×N DFT와 M×N 역 DFT의 관계는 다음과 같이 표현됩니다.

F(p,q)=m=0M1n=0N1f(m,n)ej2πpm/Mej2πqn/N   p=0, 1, ..., M1q=0, 1, ..., N1

f(m,n)=1MNp=0M1q=0N1F(p,q)ej2πpm/Mej2πqn/N   m=0, 1, ..., M1 n=0, 1, ..., N1

F(p,q)f(m,n)의 DFT 계수입니다. 영주파수 계수 F(0,0)을 보통 "DC 성분"이라고 합니다. DC는 직류를 나타내는 전기공학 용어입니다. (MATLAB®의 행렬 인덱스는 항상 0이 아닌 1에서 시작하므로 행렬 요소 f(1,1)F(1,1)은 각각 수학적 양 f(0,0)F(0,0)에 대응됩니다.)

MATLAB 함수 fft, fft2, fftn은 각각 1차원 DFT, 2차원 DFT, N차원 DFT를 계산하는 고속 푸리에 변환 알고리즘을 구현합니다. 함수 ifft, ifft2ifftn은 역 DFT를 계산합니다.

푸리에 변환과의 관계

DFT 계수 F(p,q)는 푸리에 변환 F(ω12)의 샘플입니다.

F(p,q)=F(ω1,ω2)|ω1=2πp/Mω2=2πq/Np=0,1,...,M1q=0,1,...,N1

이산 푸리에 변환 시각화하기

  1. 푸리에 변환 정의의 예제에서 사용한 함수 f(m,n)과 비슷한 행렬 f를 생성합니다. f(m,n)은 사각형 영역 안에서는 값이 1이고 영역 밖에서는 0이었다는 사실을 기억하십시오. 이진 영상을 사용하여 f(m,n)을 표현합니다.

    f = zeros(30,30);
    f(5:24,13:17) = 1;
    imshow(f,"InitialMagnification","fit")

    Binary image representation of f(m,n)

  2. 다음 명령을 사용하여 f의 30×30 DFT를 계산하고 시각화합니다.

    F = fft2(f);
    F2 = log(abs(F));
    imshow(F2,[-1 5],"InitialMagnification","fit");
    colormap(jet); colorbar

    채우기 없이 계산한 이산 푸리에 변환

    2-D plot of the 30-by-30 discrete Fourier transform of the binary rectangular function

    이 플롯은 푸리에 변환 시각화하기에서 제시된 푸리에 변환과 다릅니다. 첫째, 푸리에 변환이 보다 성기게 샘플링되었습니다. 둘째, 영주파수 계수가 기존의 중앙 위치가 아닌 왼쪽 위 코너에 표시되었습니다.

  3. 푸리에 변환의 보다 미세한 샘플링을 얻으려면 DFT를 계산할 때 f에 0 채우기를 추가하십시오. 다음 명령을 사용하여 0 채우기와 DFT 계산을 하나의 단계로 수행할 수 있습니다.

    F = fft2(f,256,256);

    위 명령은 DFT를 계산하기 전에 f가 256×256이 되도록 0 채우기를 합니다.

    imshow(log(abs(F)),[-1 5]); colormap(jet); colorbar

    채우기를 사용하여 계산한 이산 푸리에 변환

    2-D plot of the 30-by-30 discrete Fourier transform of the binary rectangular function with zero padding. The transform with zero padding has a finer frequency resolution.

  4. 여기서도 영주파수 계수가 여전히 중앙이 아닌 왼쪽 위 코너에 표시되었습니다. 이 문제는 영주파수 계수가 중앙에 오도록 F의 사분면을 바꾸는 함수 fftshift를 사용하여 해결할 수 있습니다.

    F = fft2(f,256,256);F2 = fftshift(F);
    imshow(log(abs(F2)),[-1 5]); colormap(jet); colorbar

    결과로 생성되는 플롯은 푸리에 변환 시각화하기에 표시된 플롯과 동일합니다.

푸리에 변환의 응용

이 섹션에서는 푸리에 변환의 다양한 영상 처리 관련 응용 중 몇 가지를 보여줍니다.

선형 필터의 주파수 응답

선형 필터의 임펄스 응답의 푸리에 변환은 필터의 주파수 응답을 제공합니다. 함수 freqz2는 필터의 주파수 응답을 계산하고 표시합니다. 가우스 컨벌루션 커널의 주파수 응답은 이 필터가 저주파를 통과시키고 고주파를 감쇠한다는 사실을 보여줍니다.

h = fspecial("gaussian");
freqz2(h)

가우스 필터의 주파수 응답

Mesh plot of the magnitude of the frequency response of a Gaussian filter.

선형 필터링, 필터 설계 및 주파수 응답에 대한 자세한 내용은 Design Linear Filters in the Frequency Domain 항목을 참조하십시오.

푸리에 변환을 사용하여 고속 컨벌루션 수행하기

이 예제에서는 푸리에 변환을 사용하여 두 행렬의 고속 컨벌루션을 수행하는 방법을 보여줍니다. 푸리에 변환의 주요 속성 중 하나는 두 푸리에 변환의 곱셈이 관련 공간 함수의 컨벌루션에 대응된다는 것입니다. 이 속성은 고속 푸리에 변환과 함께 고속 컨벌루션 알고리즘의 기본을 형성합니다.

참고: FFT 기반 컨벌루션 방법은 입력값이 큰 경우에 가장 자주 사용됩니다. 입력값이 작은 경우에는 일반적으로 imfilter 함수를 사용하는 것이 더 빠릅니다.

2개의 간단한 행렬 AB를 만듭니다. A는 M×N 행렬이고 B는 P×Q 행렬입니다.

A = magic(3);
B = ones(3);

AB가 적어도 (M+P-1)×(N+Q-1)이 되도록 0 채우기를 적용합니다. (ABfft2가 가장 빠르게 수행되는 크기인 2의 거듭제곱에 해당하는 크기로 0 채우기가 적용되는 경우가 많습니다.) 이 예제에서는 행렬이 8×8이 되도록 채우기를 적용합니다.

A(8,8) = 0;
B(8,8) = 0;

fft2 함수를 사용하여 AB의 2차원 DFT를 계산합니다. 두 DFT를 곱하고, ifft2 함수를 사용하여 결과값의 2차원 역 DFT를 계산합니다.

C = ifft2(fft2(A).*fft2(B));

결과값에서 0이 아닌 부분을 추출하고, 반올림 오차로 인해 생성된 허수부를 제거합니다.

C = C(1:5,1:5);
C = real(C)
C = 5×5

    8.0000    9.0000   15.0000    7.0000    6.0000
   11.0000   17.0000   30.0000   19.0000   13.0000
   15.0000   30.0000   45.0000   30.0000   15.0000
    7.0000   21.0000   30.0000   23.0000    9.0000
    4.0000   13.0000   15.0000   11.0000    2.0000

FFT 기반 상관을 수행하여 영상 특징 찾기

이 예제에서는 푸리에 변환을 사용하여 상관을 수행하는 방법을 보여줍니다. 이 연산은 컨벌루션과 밀접하게 연관되어 있습니다. 상관은 영상에서 특징을 찾는 데 사용할 수 있습니다. 이러한 맥락에서는 상관을 형판 매칭(template matching)이라고 합니다.

샘플 영상을 작업 공간으로 읽어 들입니다.

bw = imread('text.png');

영상에서 문자 "a"를 추출하여 매칭할 형판을 만듭니다. imcrop 함수의 대화형 구문을 사용하여 형판을 만들 수도 있습니다.

a = bw(32:45,88:98);

형판 영상을 180도 회전한 다음 FFT 기반 컨벌루션 기법을 사용하여 원본 영상과 형판 영상의 상관을 계산합니다. (컨벌루션 커널을 180도 회전하면 컨벌루션이 상관과 같아집니다.) 형판을 영상에 매칭하기 위해 fft2 함수와 ifft2 함수를 사용합니다. 결과로 생성되는 영상에서 밝은 피크 부분이 문자가 있는 곳에 대응됩니다.

C = real(ifft2(fft2(bw) .* fft2(rot90(a,2),256,256)));
figure
imshow(C,[]) % Scale image to appropriate display range.

Figure contains an axes object. The axes object contains an object of type image.

영상에서 형판의 위치를 보기 위해 최대 픽셀 값을 찾은 다음 이 최댓값보다 작은 임계값을 정의합니다. 이진화된 영상은 이진화된 상관 영상에서 이러한 피크를 흰색 점으로 표시합니다. (이 예제에서는 이 Figure에서 위치를 보다 쉽게 볼 수 있도록, 이진화된 영상을 팽창시켜 점의 크기를 확대합니다.)

max(C(:))
ans = 68
thresh = 60; % Use a threshold that's a little less than max.
D = C > thresh;
se = strel('disk',5);
E = imdilate(D,se);
figure
imshow(E) % Display pixels with values over the threshold.

Figure contains an axes object. The axes object contains an object of type image.

참고 항목

| |

관련 항목