이산 푸리에 변환을 이해하는 것은 쉽다. (이해만 쉽다)
웹상에서 찾을 수 있는 정보는 단순히 입력신호를 주기함수(sin, cos)의 합으로 표현할 수 있다는 것이다. (이건 쉽다)
근데 필자는 영상처리 전공이라 1D DFT는 별 관심이 없다. 그리고 1D DFT는 쉽다.
아래가 바로 2차원에서의 이산 푸리에 변환 수식이다.
(f(x,y), F(x,y) 쓰지말어라 ㅠ 헷갈린다.)
위 수식의 복소지수함수를 풀기 위해서는 오일러의 공식이 필요하다.
오일러의 공식은 $$ {e}^{ix} = \cos x + i\sin x $$ 이다.
이를 토대로 새로운 수식을 적어보면 아래와 같다.
복소수가 있으니 DFT의 결과는 2개가 나온다. (실수 부분, 허수 부분)
위 수식을 C++로 그대로 옮기면 시간 복잡도는 $$ O({n}^{4}) $$가 나오게 된다. 코드는 아래와 같다.
std::vector<cv::Mat> my_dft_n4(cv::Mat img) {
const Real pi = acosf(-1.0);
cv::Mat dft_r = cv::Mat::zeros(img.size(), Type);
cv::Mat dft_i = cv::Mat::zeros(img.size(), Type);
for (int v = 0; v < img.rows; v++) {
for (int u = 0; u < img.cols; u++) {
Real val_r = 0.0;
Real val_i = 0.0;
for (int y = 0; y < img.rows; y++) {
for (int x = 0; x < img.cols; x++) {
Real factor = 2 * pi*((static_cast<Real>(u*x) / img.cols) + (static_cast<Real>(v*y) / img.rows));
val_r += img.at<uchar>(y, x) * cos(-factor);
val_i += img.at<uchar>(y, x) * sin(-factor);
}
}
dft_r.at<Real>(v, u) = val_r;
dft_i.at<Real>(v, u) = val_i;
}
}
return { dft_r,dft_i };
}위 코드는 너무 느려서 써먹을 수 없다. 인터넷에 찾아보면 대부분 이런 식으로 짜놓는다.
실제 이미지를 한 번만 넣어보면 이게 못 써먹는다는 것을 알 텐데 ㅠㅠ
그럼, 이걸 어떻게 복잡도를 줄이냐면...
맨 처음 2D DFT의 수식을 아래와 같이 바꿀 수 있다.
지수함수는 아래의 공식이 성립한다.
그럼 지수함수의 성질에 의해 아래 공식으로 변환이 가능하다.
여기서부터 중요한 내용이다. 2D DFT는 $$ O({n}^{4}) $$ 가 아니라 $$ O({n}^{3}) $$ 으로 계산할 수 있다. row방향으로 1d dft를 구하고 column 방향으로 1d dft를 계산하여 2d dft를 계산할 수 있다.
언뜻보면 위의 수식은 풀 수 없을 것 같지만 아래 그림처럼 파란 부분은 row 부분을 계산할 때는 v와 y는 고정된 상수값이여서 sigma 밖으로 나올 수 있다.
이제 위에 노란 줄 부분이 row 방향 1d dft 이다.
아래의 수식처럼 임시공간에 중간 결과값을 계산할 수 있다.
이제 마지막 dft(u,v)를 계산하기 위한 공식은 아래와 같다.
이 역시 오일러 공식에 의해 아래와 같이 전개할 수 있다.
그런데 이 수식을 그냥 전개하면 안된다 tmp(u,v)는 복소수 도메인이기 때문이 아래의 공식을 적용한다.
그럼 최종결과는 아래와 같고
코드는 아래와 같다.
std::vector<cv::Mat> my_dft_n3(cv::Mat img) {
const Real pi = acos(-1.0);
cv::Mat tmp_r = cv::Mat::zeros(img.size(), Type); //2d temp array(real numbers)
cv::Mat tmp_i = cv::Mat::zeros(img.size(), Type); //2d temp array(imaginary numbers)
cv::Mat dft_r = cv::Mat::zeros(img.size(), Type); //2d dft array(real numbers)
cv::Mat dft_i = cv::Mat::zeros(img.size(), Type); //2d dft array(imaginary numbers)
for (int v = 0; v < img.rows; v++) {
for (int u = 0; u < img.cols; u++) {
Real val_r = 0.0;
Real val_i = 0.0;
for (int x = 0; x < img.cols; x++) {
Real factor = -2 * pi * u * x / img.cols;
val_r += img.at<uchar>(v, x) * cos(factor);
val_i += img.at<uchar>(v, x) * sin(factor);
}
tmp_r.at<Real>(v, u) = val_r;
tmp_i.at<Real>(v, u) = val_i;
}
}
for (int v = 0; v < img.rows; v++) {
for (int u = 0; u < img.cols; u++) {
Real val_r = 0.0;
Real val_i = 0.0;
for (int y = 0; y < img.rows; y++) {
Real factor = -2 * pi * v * y / img.rows;
val_r += tmp_r.at<Real>(y, u) * cos(factor) - tmp_i.at<Real>(y, u) * sin(factor);
val_i += tmp_r.at<Real>(y, u) * sin(factor) + tmp_i.at<Real>(y, u) * cos(factor);
}
dft_r.at<Real>(v, u) = val_r;
dft_i.at<Real>(v, u) = val_i;
}
}
return { dft_r,dft_i };
}Reference
Additional
https://stackoverrun.com/ko/q/2582133 (3채널 이미지에서 푸리에 변환이 불가능한 이유)
