1. 目标
本文中将会解答以下问题：
- 什么是傅里叶变换，为什么要使用它？
- 在 OpenCV 中如何实现？
- 使用以下函数：copyMakeBorder() 、merge() 、dft() 、getOptimalDFTSize() 、log() 和normalize() 。
2. 源代码
可以从这里下载samples/cpp/tutorial_code/core/discrete_fourier_transform/discrete_fourier_transform.cpp或在 OpenCV 源代码库中找到它。
这是dft() 的示例用法：
#include "opencv2/core.hpp"
#include "opencv2/imgproc.hpp"
#include "opencv2/imgcodecs.hpp"
#include "opencv2/highgui.hpp"
#include <iostream>
using namespace cv;
using namespace std;
static void help(char ** argv)
{
cout << endl
<< "This program demonstrated the use of the discrete Fourier transform (DFT). " << endl
<< "The dft of an image is taken and it's power spectrum is displayed." << endl << endl
<< "Usage:" << endl
<< argv[0] << " [image_name -- default lena.jpg]" << endl << endl;
}
int main(int argc, char ** argv)
{
help(argv);
const char* filename = argc >=2 ? argv[1] : "lena.jpg";
Mat I = imread( samples::findFile( filename ), IMREAD_GRAYSCALE);
if( I.empty()){
cout << "Error opening image" << endl;
return EXIT_FAILURE;
}
Mat padded; //expand input image to optimal size
int m = getOptimalDFTSize( I.rows );
int n = getOptimalDFTSize( I.cols ); // on the border add zero values
copyMakeBorder(I, padded, 0, m - I.rows, 0, n - I.cols, BORDER_CONSTANT, Scalar::all(0));
Mat planes[] = {Mat_<float>(padded), Mat::zeros(padded.size(), CV_32F)};
Mat complexI;
merge(planes, 2, complexI); // Add to the expanded another plane with zeros
dft(complexI, complexI); // this way the result may fit in the source matrix
// compute the magnitude and switch to logarithmic scale
// => log(1 + sqrt(Re(DFT(I))^2 + Im(DFT(I))^2))
split(complexI, planes); // planes[0] = Re(DFT(I), planes[1] = Im(DFT(I))
magnitude(planes[0], planes[1], planes[0]);// planes[0] = magnitude
Mat magI = planes[0];
magI += Scalar::all(1); // switch to logarithmic scale
log(magI, magI);
// crop the spectrum, if it has an odd number of rows or columns
magI = magI(Rect(0, 0, magI.cols & -2, magI.rows & -2));
// rearrange the quadrants of Fourier image so that the origin is at the image center
int cx = magI.cols/2;
int cy = magI.rows/2;
Mat q0(magI, Rect(0, 0, cx, cy)); // Top-Left - Create a ROI per quadrant
Mat q1(magI, Rect(cx, 0, cx, cy)); // Top-Right
Mat q2(magI, Rect(0, cy, cx, cy)); // Bottom-Left
Mat q3(magI, Rect(cx, cy, cx, cy)); // Bottom-Right
Mat tmp; // swap quadrants (Top-Left with Bottom-Right)
q0.copyTo(tmp);
q3.copyTo(q0);
tmp.copyTo(q3);
q1.copyTo(tmp); // swap quadrant (Top-Right with Bottom-Left)
q2.copyTo(q1);
tmp.copyTo(q2);
normalize(magI, magI, 0, 1, NORM_MINMAX); // Transform the matrix with float values into a
// viewable image form (float between values 0 and 1).
imshow("Input Image" , I ); // Show the result
imshow("spectrum magnitude", magI);
waitKey();
return EXIT_SUCCESS;
}
3. 代码解释
傅立叶变换将图像分解为其正弦和余弦分量。换句话说，它将图像从其空间域转换到其频域。傅里叶变换是一种方法。这个方法是任何函数都可以用无数个正弦和余弦函数之和来精确逼近。数学上二维图像傅里叶变换是：
其中用到欧拉公式如下
这里 f 是空间域中的图像值， F 是其频域中的图像值。转换的结果是复数。可以通过实部图像和复数图像或通过幅度和相位图像来显示结果。然而，在整个图像处理算法中，只有幅度图像是有趣的，因为它包含了我们需要的关于图像几何结构的所有信息。不过，如果打算以这些形式对图像进行一些修改，然后需要重新转换它，则需要保留这两种形式。
在此示例中，将展示如何计算和显示傅里叶变换的幅度图像。在数字图像的处理中数据是离散值。例如，在基本灰度图像中，值通常介于 0 和 255 之间。因此，傅里叶变换也需要是离散类型，从而导致离散傅里叶变换 ( DFT )。当需要从几何角度确定图像的结构时，都需要使用它。以下是操作步骤（假设输入是灰度图像I）：
3.1 将图像填充到最佳尺寸
DFT 的性能取决于图像大小。对于数字2、3和5的倍数的图像大小，它往往是最快的。因此，为了获得最大的性能，一种常用的方法是填充图像的边界以获得最佳尺寸大小。getOptimalDFTSize() 返回这个最佳尺寸，可以使用copyMakeBorder() 函数来填充图像的边框（附加像素初始化为零）：
Mat padded; //expand input image to optimal size
int m = getOptimalDFTSize( I.rows );
int n = getOptimalDFTSize( I.cols ); // on the border add zero values
copyMakeBorder(I, padded, 0, m - I.rows, 0, n - I.cols, BORDER_CONSTANT, Scalar::all(0));
3.2 为复值和实值分配存储位置
傅里叶变换的结果很复杂。这意味着对于每个图像值，结果是两个图像值（实部和虚部各一个）。此外，频域范围远大于其空间域范围。因此，通常至少以浮点格式存储变换结果。所以，需要将输入图像转换为浮点类型，并用另一个通道扩展它以保存复数值：
// 把原始像素值padded 转化为 float，并
Mat planes[] = {Mat_<float>(padded), Mat::zeros(padded.size(), CV_32F)};
// 分配另一个通道保存虚部
Mat complexI;
merge(planes, 2, complexI); // Add to the expanded another plane with zeros
3.3 进行离散傅里叶变换
可以进行就地计算（直接用输出覆盖输入）：
dft(complexI, complexI); // this way the result may fit in the source matrix
3.4 将实数和复数转换为幅度
复数具有实部 ( Re ) 和虚部 (Im ) 两部分。DFT 的结果是复数。DFT 的幅度为：
翻译成 OpenCV 代码：
split(complexI, planes); // planes[0] = Re(DFT(I), planes[1] = Im(DFT(I))
magnitude(planes[0], planes[1], planes[0]);// planes[0] = magnitude
Mat magI = planes[0];
3.5 切换到对数刻度
事实证明，傅立叶系数的动态范围太大而无法在屏幕上显示。结果中有非常小的和非常大的值，正常情况下不能够观察到。因为，大的值将全部变为白点，而小值将变为黑色。为了使用灰度值进行可视化，可以将线性比例转换为对数比例：
翻译成 OpenCV 代码：
magI += Scalar::all(1); // switch to logarithmic scale
log(magI, magI);
3.6 裁剪和重新排列
还记得，在第一步，扩展了图像吗？好吧，是时候裁剪掉这些填充的值了。出于可视化目的，可以重新排列结果的象限，使原点（0，0）与图像中心相对应。
// crop the spectrum, if it has an odd number of rows or columns
magI = magI(Rect(0, 0, magI.cols & -2, magI.rows & -2));
// rearrange the quadrants of Fourier image so that the origin is at the image center
int cx = magI.cols/2;
int cy = magI.rows/2;
Mat q0(magI, Rect(0, 0, cx, cy)); // Top-Left - Create a ROI per quadrant
Mat q1(magI, Rect(cx, 0, cx, cy)); // Top-Right
Mat q2(magI, Rect(0, cy, cx, cy)); // Bottom-Left
Mat q3(magI, Rect(cx, cy, cx, cy)); // Bottom-Right
Mat tmp; // swap quadrants (Top-Left with Bottom-Right)
q0.copyTo(tmp);
q3.copyTo(q0);
tmp.copyTo(q3);
q1.copyTo(tmp); // swap quadrant (Top-Right with Bottom-Left)
q2.copyTo(q1);
tmp.copyTo(q2);
3.7 标准化
出于可视化目的再次执行此操作。现在有了幅度，但是这仍然超出了图像显示范围从零到一。我们使用cv::normalize()函数将值标准化到[0,1]范围。
normalize(magI, magI, 0, 1, NORM_MINMAX); // Transform the matrix with float values into a
// viewable image form (float between values 0 and 1).
4. 结果
一个应用想法是确定图像中存在的几何方向。例如，找出文本是否是水平的？查看一些文本，会注意到文本行也形成水平线，字母形成垂直线。在傅里叶变换的情况下，也可以看到文本片段的这两个主要组成部分。这里使用水平文本图像的和旋转文本图像来演示变换结果。
横向文本，原图如下
横向文本变换结果为：
旋转文本，原图如下:
旋转文本变换结果：
可以看到频域中最有影响的分量（幅度图像上最亮的点）跟随图像上物体的几何旋转。由此可以计算偏移量并执行图像旋转以纠正最终的校准。
参考文章:
