OpenCV离散傅里叶变换

目标

我们会为以下问题寻求答案：

• 什么是傅里叶变换，为什么要使用？
• 在OpenCV中怎么做？
• 函数的使用如：cv :: copyMakeBorder（），cv :: merge（），cv :: dft（），cv :: getOptimalDFTSize（），cv :: log（）和cv :: normalize（）。

源代码

您可以从这里下载，也可以在samples/cpp/tutorial_code/core/discrete_fourier_transform/discrete_fourier_transform.cppOpenCV源代码库中找到它。

以下是cv :: dft（）的示例用法：

 #include "opencv2/core.hpp"
    2 #include "opencv2/imgproc.hpp"
    3 #include "opencv2/imgcodecs.hpp"
    4 #include "opencv2/highgui.hpp"
    5 
    6 #include <iostream>
    7 
    8 using namespace cv;
    9 using namespace std;
   10 
   11 static void help(char* progName)
   12 {
   13     cout << endl
   14         <<  "This program demonstrated the use of the discrete Fourier transform (DFT). " << endl
   15         <<  "The dft of an image is taken and it's power spectrum is displayed."          << endl
   16         <<  "Usage:"                                                                      << endl
   17         << progName << " [image_name -- default ../data/lena.jpg] "               << endl << endl;
   18 }
   19 
   20 int main(int argc, char ** argv)
   21 {
   22     help(argv[0]);
   23 
   24     const char* filename = argc >=2 ? argv[1] : "../data/lena.jpg";
   25 
   26     Mat I = imread(filename, IMREAD_GRAYSCALE);
   27     if( I.empty())
   28         return -1;
   29 
   30     Mat padded;                            //expand input image to optimal size
   31     int m = getOptimalDFTSize( I.rows );
   32     int n = getOptimalDFTSize( I.cols ); // on the border add zero values
   33     copyMakeBorder(I, padded, 0, m - I.rows, 0, n - I.cols, BORDER_CONSTANT, Scalar::all(0));
   34 
   35     Mat planes[] = {Mat_<float>(padded), Mat::zeros(padded.size(), CV_32F)};
   36     Mat complexI;
   37     merge(planes, 2, complexI);         // Add to the expanded another plane with zeros
   38 
   39     dft(complexI, complexI);            // this way the result may fit in the source matrix
   40 
   41     // compute the magnitude and switch to logarithmic scale
   42     // => log(1 + sqrt(Re(DFT(I))^2 + Im(DFT(I))^2))
   43     split(complexI, planes);                   // planes[0] = Re(DFT(I), planes[1] = Im(DFT(I))
   44     magnitude(planes[0], planes[1], planes[0]);// planes[0] = magnitude
   45     Mat magI = planes[0];
   46 
   47     magI += Scalar::all(1);                    // switch to logarithmic scale
   48     log(magI, magI);
   49 
   50     // crop the spectrum, if it has an odd number of rows or columns
   51     magI = magI(Rect(0, 0, magI.cols & -2, magI.rows & -2));
   52 
   53     // rearrange the quadrants of Fourier image  so that the origin is at the image center
   54     int cx = magI.cols/2;
   55     int cy = magI.rows/2;
   56 
   57     Mat q0(magI, Rect(0, 0, cx, cy));   // Top-Left - Create a ROI per quadrant
   58     Mat q1(magI, Rect(cx, 0, cx, cy));  // Top-Right
   59     Mat q2(magI, Rect(0, cy, cx, cy));  // Bottom-Left
   60     Mat q3(magI, Rect(cx, cy, cx, cy)); // Bottom-Right
   61 
   62     Mat tmp;                           // swap quadrants (Top-Left with Bottom-Right)
   63     q0.copyTo(tmp);
   64     q3.copyTo(q0);
   65     tmp.copyTo(q3);
   66 
   67     q1.copyTo(tmp);                    // swap quadrant (Top-Right with Bottom-Left)
   68     q2.copyTo(q1);
   69     tmp.copyTo(q2);
   70 
   71     normalize(magI, magI, 0, 1, NORM_MINMAX); // Transform the matrix with float values into a
   72                                             // viewable image form (float between values 0 and 1).
   73 
   74     imshow("Input Image"       , I   );    // Show the result
   75     imshow("spectrum magnitude", magI);
   76     waitKey();
   77 
   78     return 0;
   79 }

说明

傅立叶变换将图像分解为其窦道和余弦分量。换句话说，它会将图像从其空间域转换到其频域。这个想法是任何函数可以用无限窦和余弦函数的和来精确地近似。傅立叶变换是一种如何做到这一点的方法。数学上二维图像傅立叶变换是：

这里f是其空间域中的图像值，在其频域中是F。转换的结果是复数。可以通过实际图像和复杂图像或通过幅度和相位图像来显示这一点。然而，在整个图像处理算法中，只有幅面图像很有趣，因为它包含了我们所需要的关于图像几何结构的所有信息。然而，如果您打算以这些形式对图像进行一些修改，那么您需要重新转换它，您将需要保留这两种形式。

在本例中，我将展示如何计算和显示傅立叶变换的幅度图像。在数字图像是离散的情况下。这意味着它们可以占用给定域值的值。例如，在基本的灰度图像中，值通常在0和255之间。因此，傅里叶变换也需要是离散傅里叶变换，导致离散傅里叶变换（DFT）。当您需要从几何角度确定图像的结构时，您将需要使用它。以下是要遵循的步骤（在灰度输入图像I的情况下）：

• 将图像展开至最佳尺寸。DFT的性能取决于图像大小。对于数字二，三和五的倍数，图像尺寸趋向于最快。因此，为了获得最大的性能，通常最好将边框值填充到图像以获得具有这种特征的大小。该品种:: getOptimalDFTSize（）返回这个最佳规模，我们可以使用CV :: copyMakeBorder（）函数来扩大图像的边界：

Mat padded;                            //expand input image to optimal size
int m = getOptimalDFTSize( I.rows );
int n = getOptimalDFTSize( I.cols ); // on the border add zero pixels
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

• 将真实和复杂的值转化为大小。复数具有真实（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];

• 切换到对数刻度。原来，傅里叶系数的动态范围太大，无法显示在屏幕上。我们有一些小而高的变化值，我们不能这样观察。因此，高价值将全部作为白点，而小的则为黑色。为了将灰度值用于可视化，我们可以将我们的线性比例变换为对数：
  

翻译为OpenCV代码：

magI += Scalar::all(1);                    // switch to logarithmic scale
log(magI, magI);

• 作物和重新排列。记住，在第一步，我们扩大了形象？那么现在是摒弃新推出的价值观的时候了。为了可视化的目的，我们还可以重新排列结果的象限，使原点（零，零）对应于图像中心。

magI = magI(Rect(0, 0, magI.cols & -2, magI.rows & -2));
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);

• 规范化。这再次为可视化目的而完成。我们现在有这样的大小，但是这仍然是我们的图像显示范围从零到一。我们使用cv :: normalize（）函数将值归一化到此范围。

normalize(magI, magI, 0, 1, NORM_MINMAX); // Transform the matrix with float values into a
                                          // viewable image form (float between values 0 and 1).

结果

应用程序的想法将是确定图像中存在的几何取向。例如，让我们来看一下文本是否是横向的？看一些文字你会注意到，文本行的形式也是水平线，字母形成垂直线。在傅里叶变换的情况下，也可以看到文本片段的这两个主要组成部分。让我们使用这个水平和这个旋转的图像关于一个文本。

在水平文本的情况下：

在旋转文本的情况下：

您可以看到频域中最有影响力的组件（幅度图像上最亮点）遵循图像上对象的几何旋转。由此，我们可以计算偏移量并执行图像旋转以纠正最终的错误对准。