在OpenCV环境下写的灰度图像二维傅里叶换、幅值计算、频谱平移和将数值归一化到0到255区间的四个函数

图像处理开发需求、图像处理接私活挣零花钱，请加微信/QQ 2487872782
图像处理开发资料、图像处理技术交流请加QQ群，群号 271891601

这篇博文中包含以下四个自己根据相关数学原理写的函数：

1 灰度图像的二维傅里叶变换，函数名cv_gray_fft2()和DFT()，cv_gray_fft2()用C语言写成，DFT()用c++写成；
2 幅值计算，函数名cv_abs()
3 频谱平移, 函数名cv_gray_fft2shift()

4 将数值归一化到0到255区间，函数名cv_range_0to255()

傅里叶变换的意义这里就不多说了，大家在高等数学、信号与系统中应该都学过，不清楚的可以去查阅相关资料！这里我简单说两点！

①图像的频率是表征图像中灰度变化剧烈程序的指标；

②图像低频分量表示背景和缓慢变化区域，高频分量表示图像的边缘、细节及相关噪声。

频谱平移的作用见我的另一篇博文，链接如下：

https://blog.csdn.net/wenhao_ir/article/details/51689960

源码如下：

源码中用到的图片coins.png下载链接：https://pan.baidu.com/s/1slilbPF

#include <opencv2/opencv.hpp>  
#include <opencv2/legacy/compat.hpp> 
#include <fstream>
using namespace std;  
#pragma comment(linker, "/subsystem:\"windows\" /entry:\"mainCRTStartup\"")  

void cv_gray_fft2(IplImage *src, IplImage *dst)  //注意:此函数仅适用灰度图像,src因为是灰度图像,所以通道数要求为1,dst的属性为:IPL_DEPTH_64F,2（两通道）
{
IplImage *image_Re = 0, *image_Im = 0, *Fourier = 0;
image_Re = cvCreateImage(cvGetSize(src), IPL_DEPTH_64F, 1);  //实部
image_Im = cvCreateImage(cvGetSize(src), IPL_DEPTH_64F, 1);  //虚部
Fourier = cvCreateImage(cvGetSize(src), IPL_DEPTH_64F, 2);//2 channels,分别存储image_Re、image_Im

// Real part conversion from u8 to 64f (double)
cvConvertScale(src, image_Re);

// Imaginary part (zeros)
cvZero(image_Im);

// Join real and imaginary parts and stock them in Fourier image
cvMerge(image_Re, image_Im, 0, 0, Fourier);

// Application of the forward Fourier transform
cvDFT(Fourier, dst, CV_DXT_FORWARD);
cvReleaseImage(&image_Re);
cvReleaseImage(&image_Im);
cvReleaseImage(&Fourier);
}

void cv_abs(IplImage *src,IplImage *fft2_out, IplImage *dst)//src就是灰度图像原始数据,fft2_out要求是cv_gray_fft2的dst数据类型,即"64F, 2";dst为"64F, 1"
{
IplImage *image_Re = 0, *image_Im = 0;
image_Re = cvCreateImage(cvGetSize(src), IPL_DEPTH_64F, 1);  //实部
image_Im = cvCreateImage(cvGetSize(src), IPL_DEPTH_64F, 1);  //虚部
cvSplit(fft2_out,image_Re,image_Im,0,0);  
cvPow(image_Re,image_Re,2);                 
    cvPow(image_Im,image_Im,2); 
cvAdd(image_Re,image_Im,image_Re,NULL);
cvPow(image_Re,dst,0.5);

cvReleaseImage(&image_Re);
cvReleaseImage(&image_Im);
}

void cv_gray_fft2shift(IplImage *src)//参数是"64F, 2"
{
int nRow, nCol, i, j, cy, cx;  
double tmp13, tmp24;

//Rearrange the quadrants of Fourier image so that the origin is at the image center  
    nRow = src->height;  
    nCol = src->width;  
    cy = nRow/2; // image center  
    cx = nCol/2;  
    //CV_IMAGE_ELEM为OpenCV定义的宏，用来读取图像的像素值，这一部分就是进行中心变换  
    for( j = 0; j < cy; j++ ){  
        for( i = 0; i < cx; i++ ){  
            //中心化，将整体份成四块进行对角交换  
            tmp13 = CV_IMAGE_ELEM( src, double, j, i);  
            CV_IMAGE_ELEM( src, double, j, i) = CV_IMAGE_ELEM(  
                src, double, j+cy, i+cx);  
            CV_IMAGE_ELEM( src, double, j+cy, i+cx) = tmp13;  
  
            tmp24 = CV_IMAGE_ELEM( src, double, j, i+cx);  
            CV_IMAGE_ELEM( src, double, j, i+cx) =  
                CV_IMAGE_ELEM( src, double, j+cy, i);  
            CV_IMAGE_ELEM( src, double, j+cy, i) = tmp24;  
        }  
    }  
}

void cv_range_0to255(IplImage *src, IplImage *dst)//两个参数都是"64F, 2"
{
double m,M; 
double scale;  
double shift;  

cvMinMaxLoc(src,&m,&M,NULL,NULL); 
scale = 255/(M - m);  
shift = -m * scale;  
cvConvertScale(src,dst,scale,shift);  
}

int main()
{
int i=0;//循环变量
//从文件中加载原图  
IplImage *pSrcImage = cvLoadImage("coins.png", CV_LOAD_IMAGE_UNCHANGED);  

//创建输出的图像
IplImage *pOutImage = cvCreateImage(cvGetSize(pSrcImage), IPL_DEPTH_64F,1);

//创建中间结果图像
IplImage *cv_gray_fft2_out = cvCreateImage(cvGetSize(pSrcImage), IPL_DEPTH_64F,2);//注意这里是二通道!
IplImage *cv_abs_out = cvCreateImage(cvGetSize(pSrcImage),  IPL_DEPTH_64F, 1);

cv_gray_fft2(pSrcImage,cv_gray_fft2_out);

double watch_cv_gray_fft2_out_Re[100];//利用调试时的局部变量窗口观察cvDFTOut实部的前100个值,看与MATLAB运行的结果否相同
double watch_cv_gray_fft2_out_Im[100];//利用调试时的局部变量窗口观察cvDFTOut虚部的前100个值,看与MATLAB运行的结果否相同
for(i=0;i<10;i++)//一般情况下观察10个值就够了
{
watch_cv_gray_fft2_out_Re[i]=cvGet2D(cv_gray_fft2_out,0,i).val[0];  
}
for(i=0;i<10;i++)//一般情况下观察10个值就够了
{
watch_cv_gray_fft2_out_Im[i]=cvGet2D(cv_gray_fft2_out,0,i).val[1];  
}

cv_abs(pSrcImage,cv_gray_fft2_out,cv_abs_out);
double watch_cv_abs_out[100];//利用调试时的局部变量窗口观察cv_abs_out的前100个值,看与MATLAB运行的结果否相同
for(i=0;i<10;i++)//一般情况下观察10个值就够了
{
watch_cv_abs_out[i]=cvGet2D(cv_abs_out,0,i).val[0];  
}

cv_gray_fft2shift(cv_abs_out);
double watch_cv_gray_fft2shift_out[100];//利用调试时的局部变量窗口观察cv_gray_fft2shift_out的前100个值,看与MATLAB运行的结果否相同
for(i=0;i<10;i++)//一般情况下观察10个值就够了
{
watch_cv_gray_fft2shift_out[i]=cvGet2D(cv_abs_out,0,i).val[0];  
}

cv_range_0to255(cv_abs_out,pOutImage);
double watch_cv_range_0to255_out[100];//利用调试时的局部变量窗口观察cv_range_0to255_out的前100个值,看与MATLAB运行的结果否相同
for(i=0;i<10;i++)//一般情况下观察10个值就够了
{
watch_cv_range_0to255_out[i]=cvGet2D(pOutImage,0,i).val[0];  
}
return 0;
}

上面四个函数与MATLAB运行的结果一致,MATLAB程序如下：

clear all;
clc;
I=imread('coins.png');
F=fft2(I);
MatlabFFt2Out=F;
MatlabAbsOut=abs(F);
MatlabFftshiftOut=fftshift(MatlabAbsOut);
T=MatlabFftshiftOut;
max_T=max(max(T));
min_T=min(min(T));
shift_T=-min_T*255/(max_T-min_T);
MatlabRange0to255=T.*255/(max_T-min_T)+shift_T;

运行结果截图如下：

cv_gray_fft2

cv_abs

cv_gray_fft2shift

cv_range_0to255

下面再给出OpenCV下的C++代码实现，代码中的DFT函数已经包含了上述函数的功能：

//OpenCV版本2.4.9    
//交流QQ2487872782   
  
# include <opencv2/core/core.hpp>  
# include <opencv2/highgui/highgui.hpp>  
# include <opencv2/imgproc/imgproc.hpp>  
#include <iostream>  
using namespace cv;  
using namespace std;  
cv::Mat DFT(cv::Mat srcImage)  
{  
       cv::Mat srcGray;  
       cvtColor(srcImage,srcGray,CV_RGB2GRAY);  
       // 将输入图像延扩到最佳的尺寸(快速傅里叶变换是基于图像尺寸是2、3或5的倍数完成的,因此对于输入源图像,首先应将其转换成DFTsize，OpenCV中提供了函数getOptimalDFTSize来实现尺寸转换)  
       int nRows = getOptimalDFTSize(srcGray.rows);  
       int nCols = getOptimalDFTSize(srcGray.cols);  
       cv::Mat resultImage;  
       // 把灰度图像放在左上角,在右边和下边扩展图像,  
       // 添加的像素初始化为0  
       copyMakeBorder(srcGray, resultImage, 0,  
            nRows - srcGray.rows,   
            0, nCols - srcGray.cols,   
            BORDER_CONSTANT, Scalar::all(0));  
       // 为傅立叶变换的结果(实部和虚部)分配存储空间  
       cv::Mat planes[] = { cv::Mat_<float>(resultImage),   
            cv::Mat::zeros(resultImage.size(), CV_32F)};//这里实际上是建立MAT数组,数组有两个成员:  
														//第一个就是sizeConvMat这个对象(只是数据类型转换成了float类型)  
                                                        //第二个是全0的类型为32F的对象  
       Mat completeI;  
       // 为延扩后的图像增添一个初始化为0的通道  
       merge(planes,2,completeI); //把groupMats的第0和第1个对象合并到mergeMat,通过这个操作mergeMat是双通道的数据阵列了  
       // 进行离散傅立叶变换  
       dft(completeI,completeI);  
       // 将复数转换为幅度  
       split(completeI,planes);  
       magnitude(planes[0],planes[1],planes[0]);//0中存的是实部值,1中存的是虚部值  
       cv::Mat dftResultImage = planes[0];   
       // 对数尺度(logarithmic scale)缩放  
       dftResultImage += 1;//阵列加1作对数变换,以扩大频域动态显示范围  
       log(dftResultImage,dftResultImage);//作对数变换  
       // 剪切和重分布幅度图象限  
       dftResultImage= dftResultImage(Rect(0,  
            0,srcGray.cols,srcGray.rows));  
       // 归一化图像  
       normalize(dftResultImage,dftResultImage,  
            0,1,CV_MINMAX);  
       int cx = dftResultImage.cols/2;  
       int cy = dftResultImage.rows/2;   
       Mat tmp;  
       // Top-Left - 为每一个象限创建ROI  
       Mat q0(dftResultImage,Rect(0,0,cx,cy));  
       // Top-Right  
       Mat q1(dftResultImage,Rect(cx,0,cx,cy));  
       // Bottom-Left  
       Mat q2(dftResultImage,Rect(0,cy,cx,cy));  
       // Bottom-Right  
       Mat q3(dftResultImage,Rect(cx,cy,cx,cy));   
       // 交换象限 (Top-Left with Bottom-Right)      
       q0.copyTo(tmp);  
       q3.copyTo(q0);  
       tmp.copyTo(q3);  
       // 交换象限 (Top-Right with Bottom-Left)  
       q1.copyTo(tmp);  
       q2.copyTo(q1);  
       tmp.copyTo(q2);  
       return dftResultImage;  
}  
int main()  
{  
      cv::Mat srcImage = imread("coins.png");  
      if(srcImage.empty())  
      return-1;  
      imshow("srcImage", srcImage);  
      cv::Mat resultImage = DFT(srcImage);  
      imshow("resultImage", resultImage);  
      cv::waitKey(0);  
      return 0;  
}  

运行结果如下图所示：

图像处理开发需求、图像处理接私活挣零花钱，请加微信/QQ 2487872782
图像处理开发资料、图像处理技术交流请加QQ群，群号 271891601