关于 傅里叶变换,讲的太多了,这里我就不再啰嗦一遍了,原理的东西大家可以搜一下,推荐一篇文章:
这篇文章写得很不错了,从 频域 到傅里叶级数 讲的都比较细,而且都有配图,真看不懂的话就 你掐死他吧!
(PS:冤冤相报何时了,我在一旁看热闹,横批:不嫌事大)
看下效果(居然暴漏了自己的迅雷,哈哈):
参考代码:
/* linolzhang 2013.11
基于OpenCV的傅里叶变换及逆变换
*/
#include <iostream>
#include "opencv2/highgui/highgui.hpp"
#pragma comment(lib,"opencv_core2410.lib")
#pragma comment(lib,"opencv_highgui2410.lib")
//傅里叶正变换
void fft2(IplImage *src, IplImage *dst)
{
// 实部、虚部
IplImage *image_Re = cvCreateImage(cvGetSize(src), IPL_DEPTH_64F, 1);
IplImage *image_Im = cvCreateImage(cvGetSize(src), IPL_DEPTH_64F, 1);
IplImage *Fourier = cvCreateImage(cvGetSize(src), IPL_DEPTH_64F, 2);
cvConvertScale(src, image_Re); // 实部变换 u8 -> 64f
cvZero(image_Im); // 虚部变换 -> 0
// Merge
cvMerge(image_Re, image_Im, 0, 0, Fourier);
// DFT计算傅里叶变换
cvDFT(Fourier, dst, CV_DXT_FORWARD);
cvReleaseImage(&image_Re);
cvReleaseImage(&image_Im);
cvReleaseImage(&Fourier);
}
void fft2shift(IplImage *src, IplImage *dst)
{
IplImage *image_Re = cvCreateImage(cvGetSize(src), IPL_DEPTH_64F, 1);
IplImage *image_Im = cvCreateImage(cvGetSize(src), IPL_DEPTH_64F, 1);
cvSplit(src, image_Re, image_Im, 0, 0);
// 具体原理见冈萨雷斯《数字图像处理》p123
// Compute the magnitude of the spectrum Mag = sqrt(Re^2 + Im^2)
// 计算傅里叶谱
cvPow(image_Re, image_Re, 2.0);
cvPow(image_Im, image_Im, 2.0);
cvAdd(image_Re, image_Im, image_Re);
cvPow(image_Re, image_Re, 0.5);
// 对数变换以增强灰度级细节(这种变换使以窄带低灰度输入图像值映射
// 一宽带输出值,具体可见冈萨雷斯数字图像处理p62)
// Compute log(1 + Mag);
cvAddS(image_Re, cvScalar(1.0), image_Re); // 1 + Mag
cvLog(image_Re, image_Re); // log(1 + Mag)
// Rearrange the quadrants of Fourier image so that the origin is at the image center
int cy = src->height / 2;
int cx = src->width / 2;
// CV_IMAGE_ELEM为OpenCV定义的宏,用来读取图像的像素值,这一部分就是进行中心变换
for(int j=0; j<cy; j++)
{
for(int i=0; i<cx; i++)
{
// 中心化,将整体份成四块进行对角交换
double tmp13 = CV_IMAGE_ELEM( image_Re, double, j, i);
CV_IMAGE_ELEM( image_Re, double, j, i) = CV_IMAGE_ELEM(image_Re, double, j+cy, i+cx);
CV_IMAGE_ELEM( image_Re, double, j+cy, i+cx) = tmp13;
double tmp24 = CV_IMAGE_ELEM( image_Re, double, j, i+cx);
CV_IMAGE_ELEM( image_Re, double, j, i+cx) = CV_IMAGE_ELEM( image_Re, double, j+cy, i);
CV_IMAGE_ELEM( image_Re, double, j+cy, i) = tmp24;
}
}
// 归一化处理将矩阵的元素值归一为[0,255]
// [(f(x,y)-minVal)/(maxVal-minVal)]*255
double minVal = 0, maxVal = 0;
// Localize minimum and maximum values
cvMinMaxLoc( image_Re, &minVal, &maxVal );
// Normalize image (0 - 255) to be observed as an u8 image
double scale = 255/(maxVal - minVal);
double shift = -minVal * scale;
cvConvertScale(image_Re, dst, scale, shift);
cvReleaseImage(&image_Re);
cvReleaseImage(&image_Im);
}
int main(int argc, char** argv)
{
IplImage *src = NULL;
if(argc != 2)
src = cvLoadImage("1.jpg", 0); // 加载源图像,转为单通道
else
src = cvLoadImage(argv[1], 0);
IplImage *Fourier = cvCreateImage(cvGetSize(src),IPL_DEPTH_64F,2); // 傅里叶系数
IplImage *dst = cvCreateImage(cvGetSize(src),IPL_DEPTH_64F,2);
IplImage *ImageRe = cvCreateImage(cvGetSize(src),IPL_DEPTH_64F,1);
IplImage *ImageIm = cvCreateImage(cvGetSize(src),IPL_DEPTH_64F,1);
IplImage *Image = cvCreateImage(cvGetSize(src),src->depth,src->nChannels);
IplImage *ImageDst = cvCreateImage(cvGetSize(src),src->depth,src->nChannels);
fft2(src,Fourier); // 傅里叶变换
fft2shift(Fourier, Image); // 中心化
cvDFT(Fourier,dst,CV_DXT_INV_SCALE); // 实现傅里叶逆变换,并对结果进行缩放
cvSplit(dst,ImageRe,ImageIm,0,0);
// 对数组每个元素平方并存储在第二个参数中
cvPow(ImageRe,ImageRe,2);
cvPow(ImageIm,ImageIm,2);
cvAdd(ImageRe,ImageIm,ImageRe,NULL);
cvPow(ImageRe,ImageRe,0.5);
double m,M;
cvMinMaxLoc(ImageRe,&m,&M,NULL,NULL);
double scale = 255/(M - m);
double shift = -m * scale;
// 将shift加在ImageRe各元素按比例缩放的结果上,存储为ImageDst
cvConvertScale(ImageRe,ImageDst,scale,shift);
cvShowImage("源图像",src);
cvShowImage("傅里叶谱",Image);
cvShowImage("傅里叶逆变换",ImageDst);
cvWaitKey(0);
// 释放资源
cvReleaseImage(&src);
cvReleaseImage(&Image);
cvReleaseImage(&ImageIm);
cvReleaseImage(&ImageRe);
cvReleaseImage(&Fourier);
cvReleaseImage(&dst);
cvReleaseImage(&ImageDst);
return 0;
}