首页 > 代码库 > OpenCV Tutorials —— Discrete Fourier Transform

OpenCV Tutorials —— Discrete Fourier Transform

The Fourier Transform will decompose an image into its sinus and cosines components. In other words, it will transform an image from its spatial domain to its frequency domain.

将图像从空域转换到频域,使其由 sin 和 cos 成分构成

The idea is that any function may be approximated exactly with the sum of infinite sinus and cosines functions.

F(k,l) = \displaystyle\sum\limits_{i=0}^{N-1}\sum\limits_{j=0}^{N-1} f(i,j)e^{-i2\pi(\frac{ki}{N}+\frac{lj}{N})}

e^{ix} = \cos{x} + i\sin {x}

The result of the transformation is complex numbers. Displaying this is possible either via a real image and a complex image or via a magnitude and a phaseimage.

实数图像&复数图像    或者     幅度谱&相位谱

 

In case of digital images are discrete. This means they may take up a value from a given domain value. For example in a basic gray scale image values usually are between zero and 255. Therefore the Fourier Transform too needs to be of a discrete type resulting in a Discrete Fourier Transform (DFT). You’ll want to use this whenever you need to determine the structure of an image from a geometrical point of view.

Expand the image to an optimal size.It tends to be the fastest for image sizes that are multiple of the numbers two, three and five. Therefore, to achieve maximal performance it is generally a good idea to pad border values to the image to get a size with such traits. ThegetOptimalDFTSize() returns this optimal size and we can use the copyMakeBorder() function to expand the borders of an image:

返回最优DFT参数     扩展图像边界

 

Make place for both the complex and the real values The result of a Fourier Transform is complex. This implies that for each image value the result is two image values (one per component). Moreover, the frequency domains range is much larger than its spatial counterpart. Therefore, we store these usually at least in a float format.

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

 

Make the Discrete Fourier Transform. It’s possible an in-place calculation (same input as output)

dft(complexI, complexI);            // this way the result may fit in the source matrix

 

Transform the real and complex values to magnitude. A complex number has a real (Re) and a complex (imaginary -Im) part. The results of a DFT are complex numbers. The magnitude of a DFT is:

M = \sqrt[2]{ {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] = magnitudeMat magI = planes[0];

 

Switch to a logarithmic scale. It turns out that the dynamic range of the Fourier coefficients is too large to be displayed on the screen. We have some small and some high changing values that we can’t observe like this. Therefore the high values will all turn out as white points, while the small ones as black. To use the gray scale values to for visualization we can transform our linear scale to a logarithmic one:

M_1 = \log{(1 + M)}

 

magI += Scalar::all(1);                    // switch to logarithmic scalelog(magI, magI);	// 对数不能为零

 

Crop and rearrange. Remember, that at the first step, we expanded the image? Well, it’s time to throw away the newly introduced values. For visualization purposes we may also rearrange the quadrants of the result, so that the origin (zero, zero) corresponds with the image center.

切掉第一步扩展的部分

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 quadrantMat q1(magI, Rect(cx, 0, cx, cy));  // Top-RightMat q2(magI, Rect(0, cy, cx, cy));  // Bottom-LeftMat q3(magI, Rect(cx, cy, cx, cy)); // Bottom-RightMat 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);

 

  1. Normalize. This is done again for visualization purposes. We now have the magnitudes, however this are still out of our image display range of zero to one. We normalize our values to this range using the normalize() function.

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

 

注意:

Mat_<float>(padded) 强制类型转换

merge —— split

magnitude  —— 计算模长

The function magnitude calculates the magnitude of 2D vectors formed from the corresponding elements of x and yarrays:

\texttt{dst} (I) =  \sqrt{\texttt{x}(I)^2 + \texttt{y}(I)^2}

 

copyTo

C++: void Mat::copyTo(OutputArray m) const
C++: void Mat::copyTo(OutputArray m, InputArray mask) const

Parameters:

  • m – Destination matrix. If it does not have a proper size or type before the operation, it is reallocated.
  • mask – Operation mask. Its non-zero elements indicate which matrix elements need to be copied.

 

QTDFBO7)W187{BH[X1]$QEK

ALL Code

#include "stdafx.h"#include "opencv2/core/core.hpp"#include "opencv2/imgproc/imgproc.hpp"#include "opencv2/highgui/highgui.hpp"#include <iostream>using namespace cv;using namespace std;static void help(char* progName){	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		<<  "Usage:"                                                                      << endl		<< progName << " [image_name -- default lena.jpg] "                       << endl << endl;}int main(int argc, char ** argv){	// help(argv[0]);	const char* filename = argc >=2 ? argv[1] : "xue.jpg";	Mat I = imread(filename, CV_LOAD_IMAGE_GRAYSCALE);	if( I.empty())		return -1;	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);		// 都是mat引用,指向图像的不同区域 ~~ 就像两个变量值互换	normalize(magI, magI, 0, 1, CV_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 0;}

OpenCV Tutorials —— Discrete Fourier Transform