傅里叶变换及低通滤波再反变换（C++&&opencv）

对灰度图像做傅里叶变换并进行低通滤波再做傅里叶逆变换。

结果及代码如下：

 如上图：mag图表示幅度谱，亮处为低频部分，暗处为高频部分。

#include <iostream>
#include <opencv.hpp>
#include <opencv2\opencv.hpp>
#include <opencv2/highgui/highgui.hpp>    
#include <opencv2/core/core.hpp>
#include <opencv2/highgui/highgui_c.h>

using namespace cv;
using namespace std;

void FourierTransform(cv::Mat& image)
{
	image.convertTo(image, CV_32F);
	vector<Mat> channels;
	split(image, channels);  //分离RGB通道
	Mat image_B = channels[0];
	//选取最适合做fft的宽和高
	int m1 = getOptimalDFTSize(image_B.rows);  
	int n1 = getOptimalDFTSize(image_B.cols);
	Mat padded;
	//填充
	copyMakeBorder(image_B, padded, 0, m1 - image_B.rows, 0, n1 - image_B.cols, BORDER_CONSTANT, Scalar::all(0));
	Mat planes[] = { Mat_<float>(padded), Mat::zeros(padded.size(), CV_32F) };
	Mat complexI;
	merge(planes, 2, complexI);  //planes[0], planes[1]是实部和虚部

	dft(complexI, complexI, DFT_SCALE | DFT_COMPLEX_OUTPUT);
	split(complexI, planes);

	//定义幅度谱和相位谱
	Mat ph, mag, idft;
	phase(planes[0], planes[1], ph);
	magnitude(planes[0], planes[1], mag);  //由实部planes[0]和虚部planes[1]得到幅度谱mag和相位谱ph

	//重新排列傅里叶图像中的象限，使得原点位于图像中心
	int cx = mag.cols / 2;
	int cy = mag.rows / 2;
	Mat q0(mag, Rect(0, 0, cx, cy));       //左上角图像划定ROI区域
	Mat q1(mag, Rect(cx, 0, cx, cy));      //右上角图像
	Mat q2(mag, Rect(0, cy, cx, cy));      //左下角图像
	Mat q3(mag, Rect(cx, cy, cx, cy));     //右下角图像

	//变换左上角和右下角象限
	Mat tmp;
	q0.copyTo(tmp);
	q3.copyTo(q0);
	tmp.copyTo(q3);

	//变换右上角和左下角象限
	q1.copyTo(tmp);
	q2.copyTo(q1);
	tmp.copyTo(q2);

	imshow("mag", mag);

	//低通滤波
	for (int i = 0; i < mag.cols;i++){
		for (int j = 0; j < mag.rows; j++){
			if (abs(i - mag.cols / 2) > mag.cols / 10 || abs(j - mag.rows / 2) > mag.rows / 10)
				mag.at<float>(j, i) = 0;
		}
	}

	imshow("mag2", mag);
	q0.copyTo(tmp);
	q3.copyTo(q0);
	tmp.copyTo(q3);

	//变换右上角和左下角象限
	q1.copyTo(tmp);
	q2.copyTo(q1);
	tmp.copyTo(q2);

	//傅里叶逆变换
	polarToCart(mag, ph, planes[0], planes[1]);  //由幅度谱mag和相位谱ph恢复实部planes[0]和虚部planes[1]
	merge(planes, 2, idft);
	dft(idft, idft, DFT_INVERSE | DFT_REAL_OUTPUT);
	image_B = idft(Rect(0, 0, image.cols & -2, image.rows & -2));
	image_B.copyTo(channels[0]);
	
	merge(channels, image);
	image.convertTo(image, CV_8U);
}


void main()
{
	Mat img = imread("E://dx.jpg");
	cvtColor(img, img, COLOR_BGR2GRAY);
	imshow("src", img);
	FourierTransform(img);
	imshow("DFT img", img);
	waitKey();
	system("pause");
	return ;
}

附鼠标响应的傅里叶变换并滤波代码

 

#include <iostream>
#include <opencv.hpp>
#include <opencv2\opencv.hpp>
#include <opencv2/highgui/highgui.hpp>    
#include <opencv2/core/core.hpp>
#include <opencv2/highgui/highgui_c.h>

using namespace cv;
using namespace std;

struct fft 
{
	Mat img;
	vector<Point2i> filt_points;
	Point2i thisPoint;
	int filt_step;
	Mat src_mag, ph, dst;
	Mat planes[2];
};

void IDFTtransform(Mat& image,Mat& src_mag, Mat& ph, Mat planes[])
{
	imshow("src_mag", src_mag);

	vector<Mat> channels;
	split(image, channels);  //分离RGB通道
	Mat image_B = channels[0];

	int cx = src_mag.cols / 2;
	int cy = src_mag.rows / 2;
	Mat q0(src_mag, Rect(0, 0, cx, cy));       //左上角图像划定ROI区域
	Mat q1(src_mag, Rect(cx, 0, cx, cy));      //右上角图像
	Mat q2(src_mag, Rect(0, cy, cx, cy));      //左下角图像
	Mat q3(src_mag, Rect(cx, cy, cx, cy));     //右下角图像

	//变换左上角和右下角象限
	Mat tmp,idft;
	q0.copyTo(tmp);
	q3.copyTo(q0);
	tmp.copyTo(q3);

	//变换右上角和左下角象限
	q1.copyTo(tmp);
	q2.copyTo(q1);
	tmp.copyTo(q2);

	polarToCart(src_mag, ph, planes[0], planes[1]);  //由幅度谱mag和相位谱ph恢复实部planes[0]和虚部planes[1]
	merge(planes, 2, idft);
	dft(idft, idft, DFT_INVERSE | DFT_REAL_OUTPUT);
	image_B = idft(Rect(0, 0, image.cols & -2, image.rows & -2));
	image_B.copyTo(channels[0]);

	merge(channels, image);
	image.convertTo(image, CV_8U);

	q0.copyTo(tmp);
	q3.copyTo(q0);
	tmp.copyTo(q3);

	//变换右上角和左下角象限
	q1.copyTo(tmp);
	q2.copyTo(q1);
	tmp.copyTo(q2);
	imshow("dst", image);
}

//DFTtransform实现
vector<Mat> DFTtransform(Mat& image, Mat& src_mag, Mat& ph, Mat& dst)
{
	image.convertTo(image, CV_32F);
	vector<Mat> channels;
	split(image, channels);  //分离RGB通道
	Mat image_B = channels[0];
	//expand input image to optimal size
	int m1 = getOptimalDFTSize(image_B.rows);  //选取最适合做fft的宽和高
	int n1 = getOptimalDFTSize(image_B.cols);
	Mat padded;
	//填充0
	copyMakeBorder(image_B, padded, 0, m1 - image_B.rows, 0, n1 - image_B.cols, BORDER_CONSTANT, Scalar::all(0));
	Mat planes[] = { Mat_<float>(padded), Mat::zeros(padded.size(), CV_32F) };
	Mat complexI;
	merge(planes, 2, complexI);  //planes[0], planes[1]是实部和虚部

	dft(complexI, complexI, DFT_SCALE | DFT_COMPLEX_OUTPUT);
	split(complexI, planes);

	//定义幅度谱和相位谱
	Mat  mag, idft;
	phase(planes[0], planes[1], ph);
	magnitude(planes[0], planes[1], src_mag);  //由实部planes[0]和虚部planes[1]得到幅度谱mag和相位谱ph

	//重新排列傅里叶图像中的象限，使得原点位于图像中心
	int cx = src_mag.cols / 2;
	int cy = src_mag.rows / 2;
	Mat q0(src_mag, Rect(0, 0, cx, cy));       //左上角图像划定ROI区域
	Mat q1(src_mag, Rect(cx, 0, cx, cy));      //右上角图像
	Mat q2(src_mag, Rect(0, cy, cx, cy));      //左下角图像
	Mat q3(src_mag, Rect(cx, cy, cx, cy));     //右下角图像

	//变换左上角和右下角象限
	Mat tmp;
	q0.copyTo(tmp);
	q3.copyTo(q0);
	tmp.copyTo(q3);

	//变换右上角和左下角象限
	q1.copyTo(tmp);
	q2.copyTo(q1);
	tmp.copyTo(q2);

	vector<Mat> plane12;
	plane12.push_back(planes[0]);
	plane12.push_back(planes[1]);
	return (plane12);
}

void on_mouse(int event, int x, int y, int flags, void* lst)
{
	fft* temp = (fft*)lst;
	if (event == CV_EVENT_LBUTTONDOWN){
		cout << "x:" << x << " y:" << y << endl;
		temp->thisPoint.x = x;
		temp->thisPoint.y = y;
		temp->filt_points.push_back(temp->thisPoint);
		for (int i = x - 5; i < x + 6;i++){
			for (int j = y - 5; j < y + 6; j++){
				temp->src_mag.at<float>(j, i) = 0;
			}
		}
		IDFTtransform(temp->img, temp->src_mag, temp->ph, temp->planes);
		temp->filt_step++;
	}
}

void main(int event, int x, int y, int flags, void *ustc)
{
	fft lst;
	lst.filt_step=0;
	lst.img = imread("E://dx.jpg");
	cvtColor(lst.img, lst.img, COLOR_BGR2GRAY);
	imshow("src", lst.img);
	vector<Mat> planes12;

	planes12 = DFTtransform(lst.img, lst.src_mag, lst.ph, lst.dst);
	lst.planes[0] = planes12[0]; 
	lst.planes[1] = planes12[1];
	IDFTtransform(lst.img, lst.src_mag, lst.ph, lst.planes);

	Point2i thisPoint;
	setMouseCallback("src_mag", on_mouse, (void*)&lst);

	waitKey();
	system("pause");
	return ;
}