opencv mat初始化_【2】OpenCV核心模块（7）DFT离散傅里叶变换

本篇讲解图像的离散傅里叶变换DFT。通过DFT我们可以获取图像的频域信息，根据频谱能够获取图像的几何结构特性。本节利用OpenCV提供的一系列函数实现DFT，并显示了结果。最后，介绍了DFT在旋转文本矫正中的作用。用到了下面六个函数。copyMakeBorder() , merge() , dft() , getOptimalDFTSize() , log() 和 normalize() 。

• 原文网址Discrete Fourier Transform
• 本地目录D:opencvsourcesdoctutorialscorediscrete_fourier_transform
• 代码目录D:opencvsourcessamplescpptutorial_codecorediscrete_fourier_transform
• GitHub 有相应文档和OpenCV源代码
• 版本OpenCV4.1.2（版本兼容性见英文原文，部分文档适用于OpenCV2.0和3.0）
• 环境Windows、C++、VS2019 Community

---

代码

下面代码是函数dft()的用法示例：

//---图像离散傅里叶变换DFT，显示变换结果---//
//---修改CVer---------------------//
#include <opencv2opencv.hpp>
#include <iostream>

#ifdef _DEBUG
#pragma comment(lib,"opencv_world412d.lib")
#else
#pragma comment(lib,"opencv_world412.lib")
#endif

#include <iostream>

using namespace cv;
using namespace std;
		
int main(int argc, char** argv)
{
	const char* filename = "imageTextN.png"; //OpenCV提供的文本图像，D:opencvsourcessamplesdata
	//const char* filename = "imageTextR.png";//OpenCV提供的旋转文本图像
          
	Mat I = imread(samples::findFile(filename), IMREAD_GRAYSCALE);//读入图像为灰度
	if (I.empty()) 
	{
		cout << "Error opening image" << endl;
		return EXIT_FAILURE;
	}

	//优化输入图像大小，对边界补零，便于DFT计算
	Mat padded;                            
	int m = getOptimalDFTSize(I.rows);
	int n = getOptimalDFTSize(I.cols);     
	copyMakeBorder(I, padded, 0, m - I.rows, 0, n - I.cols, BORDER_CONSTANT, Scalar::all(0));//图像边界补零
	imshow("", padded); waitKey(0);

	//创建Mat数组planes，两个元素，一个是扩展后并转为float型的图像，另一个是同样大小，赋值为0的图像。
	Mat planes[] = { Mat_<float>(padded), Mat::zeros(padded.size(), CV_32F) };
	Mat complexI;
	merge(planes, 2, complexI);         // 合并Mat数组的两个元素，形成一个双通道图像complexI

	dft(complexI, complexI);            //  dft运算，结果仍然存放在complexI

	// 计算幅值并转换到对数尺度，用下面公式
	// => log(1 + sqrt(Re(DFT(I))^2 + Im(DFT(I))^2))
	split(complexI, planes);                   // 分离complexI的两个通道到数组 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); // 将变换结果转换到0~1之间，便于显示
	imshow("Input Image", I);    // 显示原图Show the result
	imshow("spectrum magnitude", magI);//显示结果
	waitKey();

	return EXIT_SUCCESS;
}

解释

图像经过傅里叶变换会分解为正弦和余弦成分。换句话说，会把一副图像从空间域转换到频域。任意函数都可能被无限的正弦与余弦函数的和精确逼近，傅里叶变换可以实现这种逼近。一副两维图像傅里叶变换的数学表达如下：

这里f是图像空域的值（像素值），F是频域的。变换的结果是复数。如果想观察这个示结果，既可以观察实数部分以及虚数部分的图像，也可以观察变换结果的幅度和相位图像。然而，整个算法处理过程，只有幅度图像是我们感兴趣的，因为它包含图像几何结构的所有信息。但是，如果想对变换结果进行修改，则需要用到所有变换结果。

这个例子展示了如何计算和显示傅里叶变换后的幅度图像。因为数字图像是离散的，像素的取值是离散的，因此这里的傅里叶变换是离散形式，即DFT。如果想从几何角度了解图像结构就会用到这个方法。下面是步骤（以灰度图像I为例）：

将图像扩展到最佳尺寸 Expand the image to an optimal size

DFT的性能与输入图像的尺寸有关。如果图像的尺寸是2、3、5的乘积构成，速度最快。因此，为了实现最大性能，通常对图像边界进行补充，使得能够满足最佳尺寸要求。使用函数getOptimalDFTSize()会返回这个最佳尺寸，然后使用copyMakeBorder()函数来扩展图像 (扩展部分的像素值初始化为0)：

//优化输入图像大小，对边界补零，便于DFT计算
Mat padded;                            
int m = getOptimalDFTSize(I.rows);
int n = getOptimalDFTSize(I.cols);     
copyMakeBorder(I, padded, 0, m - I.rows, 0, n - I.cols, BORDER_CONSTANT, Scalar::all(0));//图像边界补零

创建变量存放实部和虚部 Make place for both the complex and the real values

傅里叶变换的结果是复数。这意味着每幅图形的变换结果是两幅图像（复数的实部和虚部）。另外，频域取值范围远大于对应的空域部分。因此，变换结果至少用float类型的数据。所以输入图像需要转化为该类型，并且增加一个通道存放虚部。

//创建Mat数组planes，两个元素，一个是扩展后并转为float型的图像，另一个是同样大小，赋值为0的图像。
Mat planes[] = { Mat_<float>(padded), Mat::zeros(padded.size(), CV_32F) };
Mat complexI;
merge(planes, 2, complexI);         // 合并Mat数组的两个元素，形成一个双通道图像complexI

执行DFT函数 Make the Discrete Fourier Transform

这是一个in-place计算(输入输出用同一个变量):

dft(complexI, complexI);            //  dft运算，结果仍然存放在complexI

将实部和虚部转换为幅值 Transform the real and complex values to magnitude

一个复数有实部（Real - Re）和虚部 (imaginary - Im) 。DFT的结果是复数。幅值用下面公式计算：

// 计算幅值并转换到对数尺度，用下面公式
// => log(1 + sqrt(Re(DFT(I))^2 + Im(DFT(I))^2))
split(complexI, planes);                   // 分离complexI的两个通道到数组 planes[0] = Re(DFT(I), planes[1] = Im(DFT(I))
magnitude(planes[0], planes[1], planes[0]);// 数组的第一个元素是幅度图像 planes[0] = magnitude
Mat magI = planes[0];

转换到对数刻度Switch to a logarithmic scale

傅里叶系数的动态范围太大不适宜用图像显示，一些较小和快速变化的值不容易观察到。因此， 一些大的值是白点，而小的是黑点。便于观察图像，进行下面变换，显示为灰度图像：

magI += Scalar::all(1);                    // 幅值转换到对数尺度switch to logarithmic scale
log(magI, magI);

截取并重新组合图像 Crop and rearrange

将幅度图像转换为偶数。并且重新组合图像象限，让原点（0,0）位于图像中心。

// 重新排列变换结果，让原点处于图像中心。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

上面经过对数变换后，数值可能仍然不在0和1之间。这里使用cv::normalize()函数进行变换。

normalize(magI, magI, 0, 1, NORM_MINMAX); // 将变换结果转换到0~1之间，便于显示

结果 Result

有一种图像处理需求是找出图像中的几何结构方向。例如，一副图像中的文本是否水平。常见的文本是按行排列的，构成水平线。而其中的字符可能形成了垂直线。文本的这两种成分通过傅里叶变换可以观察到。下面分别是是水平和旋转后文本图像的傅里叶变换结果，（可以看到两幅输入图像进行了补零扩展，根据图像大小不同，补零情况也不同）。

正常文本及其傅里叶变换结果

旋转文本及其傅里叶变换

从上面的图中可以看到较亮的点就是起决定性的频率成分。第二幅图中，出现了一条斜线，而这条线与旋转的文本方向正好垂直。如果能提取这条线，并计算出其角度就可以对倾斜的文本进行矫正（后面图像处理部分会讲解相关的一些方法）。