l 频率域滤波的基本步骤
思想:通过滤波器函数以某种方式来修改图像 变换,然后通过取结果的反变换来获得处理后 的输出图像
低通滤波器:
使低频通过而使高频衰减的滤波 器
被低通滤波的图像比原始图像少尖锐的细节部分而 突出平滑过渡部分
对比空间域滤波的平滑处理,如均值滤波器
高通滤波器:
使高频通过而使低频衰减的滤波 器
被高通滤波的图像比原始图像少灰度级的平滑过渡 而突出边缘等细节部分
对比空间域的梯度算子、拉普拉斯算子
从滤波器处理效果的尖锐程度,可以将他们分为三种类型:理想滤波器、巴特沃斯滤波器、高斯滤波器。他们的尖锐程度也是依次递减。
再从滤波的通过范围看,这三种滤波器都有低通、高通、带通、带阻四个版本。
下面将他们的的公式贴出来,就一目了然啦。
理想低通滤波器:
理想高通滤波器:
理想带阻滤波器:
这里:
巴特沃斯低通滤波器:
巴特沃斯高通滤波器:
巴特沃斯带阻滤波器:
W是带宽
高斯低通滤波器
高斯高通滤波器
高斯带阻滤波器
按照公式可以实现任意滤波器。
频率域滤波——以高斯低通为例
//**************************************
//频率域滤波——以高斯低通为例
//**************************************
#include<opencv2/opencv.hpp>
#include<iostream>
using namespace std;
using namespace cv;
Mat gaussianlbrf(Mat scr,float sigma);//高斯低通滤波器函数
Mat freqfilt(Mat scr,Mat blur);//频率域滤波函数
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 ../data/lena.jpg] " << endl << endl;
}
int main( int argc, char *argv[])
{
help(argv[0]);
const char* filename = argc >=2 ? argv[1] : "../data/lena.jpg";
Mat input = imread(filename, IMREAD_GRAYSCALE);
if( input.empty())
return -1;
imshow("input",input);//显示原图
int w=getOptimalDFTSize(input.cols); //获取进行dtf的最优尺寸
int h=getOptimalDFTSize(input.rows); //获取进行dtf的最优尺寸
Mat padded;
copyMakeBorder(input,padded,0,h-input.rows,0,w-input.cols, BORDER_CONSTANT , Scalar::all(0)); //边界填充
padded.convertTo(padded,CV_32FC1); //将图像转换为flaot型
Mat gaussianKernel=gaussianlbrf(padded,45);//高斯低通滤波器
Mat out=freqfilt(padded,gaussianKernel);//频率域滤波
out = out(Rect(0,0,input.cols,input.rows));
imshow("结果图 sigma=40",out);
waitKey(0);
return 0;
}
//*****************高斯低通滤波器***********************
Mat gaussianlbrf(Mat scr,float sigma)
{
Mat gaussianBlur(scr.size(),CV_32FC1); //,CV_32FC1
float d0=2*sigma*sigma;//高斯函数参数,越小,频率高斯滤波器越窄,滤除高频成分越多,图像就越平滑
for(int i=0;i<scr.rows ; i++ )
{
for(int j=0; j<scr.cols ; j++ )
{
float d=pow(float(i-scr.rows/2),2)+pow(float(j-scr.cols/2),2);//分子,计算pow必须为float型
gaussianBlur.at<float>(i,j)=expf(-d/d0);//expf为以e为底求幂(必须为float型)
}
}
imshow("高斯低通滤波器",gaussianBlur);
return gaussianBlur;
}
//*****************频率域滤波*******************
Mat freqfilt(Mat scr,Mat blur)
{
//***********************DFT*******************
Mat plane[]={scr, Mat::zeros(scr.size() , CV_32FC1)}; //创建通道,存储dft后的实部与虚部(CV_32F,必须为单通道数)
Mat complexIm;
merge(plane,2,complexIm);//合并通道 (把两个矩阵合并为一个2通道的Mat类容器)
dft(complexIm,complexIm);//进行傅立叶变换,结果保存在自身
//***************中心化********************
split(complexIm,plane);//分离通道(数组分离)
int cx=plane[0].cols/2;int cy=plane[0].rows/2;//以下的操作是移动图像 (零频移到中心)
Mat part1_r(plane[0],Rect(0,0,cx,cy)); //元素坐标表示为(cx,cy)
Mat part2_r(plane[0],Rect(cx,0,cx,cy));
Mat part3_r(plane[0],Rect(0,cy,cx,cy));
Mat part4_r(plane[0],Rect(cx,cy,cx,cy));
Mat temp;
part1_r.copyTo(temp); //左上与右下交换位置(实部)
part4_r.copyTo(part1_r);
temp.copyTo(part4_r);
part2_r.copyTo(temp); //右上与左下交换位置(实部)
part3_r.copyTo(part2_r);
temp.copyTo(part3_r);
Mat part1_i(plane[1],Rect(0,0,cx,cy)); //元素坐标(cx,cy)
Mat part2_i(plane[1],Rect(cx,0,cx,cy));
Mat part3_i(plane[1],Rect(0,cy,cx,cy));
Mat part4_i(plane[1],Rect(cx,cy,cx,cy));
part1_i.copyTo(temp); //左上与右下交换位置(虚部)
part4_i.copyTo(part1_i);
temp.copyTo(part4_i);
part2_i.copyTo(temp); //右上与左下交换位置(虚部)
part3_i.copyTo(part2_i);
temp.copyTo(part3_i);
//*****************滤波器函数与DFT结果的乘积****************
Mat blur_r,blur_i,BLUR;
multiply(plane[0], blur, blur_r); //滤波(实部与滤波器模板对应元素相乘)
multiply(plane[1], blur,blur_i);//滤波(虚部与滤波器模板对应元素相乘)
Mat plane1[]={blur_r, blur_i};
merge(plane1,2,BLUR);//实部与虚部合并
//*********************得到原图频谱图***********************************
magnitude(plane[0],plane[1],plane[0]);//获取幅度图像,0通道为实部通道,1为虚部,因为二维傅立叶变换结果是复数
plane[0]+=Scalar::all(1); //傅立叶变o换后的图片不好分析,进行对数处理,结果比较好看
log(plane[0],plane[0]); // float型的灰度空间为[0,1])
normalize(plane[0],plane[0],1,0,CV_MINMAX); //归一化便于显示
imshow("原图像频谱图",plane[0]);
//******************IDFT*******************************
/*
Mat part111(BLUR,Rect(0,0,cx,cy)); //元素坐标(cx,cy)
Mat part222(BLUR,Rect(cx,0,cx,cy));
Mat part333(BLUR,Rect(0,cy,cx,cy));
Mat part444(BLUR,Rect(cx,cy,cx,cy));
part111.copyTo(temp); //左上与右下交换位置(虚部)
part444.copyTo(part111);
temp.copyTo(part444);
part222.copyTo(temp); //右上与左下交换位置
part333.copyTo(part222);
temp.copyTo(part333);
*/
idft( BLUR, BLUR); //idft结果也为复数
split(BLUR,plane);//分离通道,主要获取通道
magnitude(plane[0],plane[1],plane[0]); //求幅值(模)
normalize(plane[0],plane[0],1,0,CV_MINMAX); //归一化便于显示
return plane[0];//返回参数
}
结果:
原图:
原图像频谱图
高斯低通滤波器
结果图像:
结论
ü GLPF不能达到有相同截止频率的二阶 BLPF的平滑效果
ü GLPF没有振铃
ü 如果需要严格控制低频和高频之间截 至频率的过渡,选用BLPF,代价是可能 产生振铃