1 巴特沃斯高通滤波
clc;
clear;
close all;
%% 读取图像
Img = imread('photo.jpg');
%% 可视化原始图像
subplot(121),imshow(Img);
title('原始图像');
%% --------------------------巴特沃斯高通滤波
% 图像的数据类型由 uint8 转化为 double 类型
Imgd = im2double(Img);
% 二维傅立叶变换
fuv = fft2(Imgd);
% 对傅里叶变换后得到的频谱进行平移,将变换后的图像频谱中心从矩阵的原点移
到矩阵的中心
trans = fftshift(fuv);
[row,col] = size(trans); % 求二维傅里叶变换后图像大小
n = 2; % 巴特沃斯高通滤波器的阶数 n
d0 = 30; % 巴特沃斯高通滤波器的截止频率 D0
n1 = round(row/2);
n2 = round(col/2);
% 计算频率点(i,j)与频域中心的距离
for i = 1:row
for j = 1:col
distance = sqrt((i-n1)^2+(j-n2)^2);
if distance==0
h = 0;
else
h = 1/(1+(d0/distance)^(2*n));% 巴特沃斯高通滤波
end
trans(i,j) = h*trans(i,j);% 频域图像乘以滤波器的系数
end
end
invTrans = ifftshift(trans);
% 二维傅里叶反变换转换为时域图像
fxy = ifft2(invTrans);
% real 函数取元素的实部
trans = real(fxy);
%% ----------------------可视化滤波后的图像
--------------------------------
subplot(122),imshow(trans,[]);
title('巴特沃斯高通滤波图像');
2 巴特沃斯低通滤波
close all;
clear all;
clc;
I = imread('photo.jpg');
subplot(131),imshow(I);
title('原始图像');
I=imnoise(I,'gaussian');%加入高斯白噪声
subplot(132),imshow(I);
title('加入噪声后所得图像');
s=fftshift(fft2(im2double(I)));
[N1,N2]=size(s);%求二维傅里叶变换后图像大小
n=2; % 将巴特沃斯低通滤波器的阶数 n 设置为 2
d0=30; % 将巴特沃斯低通滤波器的截止频率 D0 设置为 30
n1=round(N1/2);
n2=round(N2/2);
for i=1:N1 %双重 for 循环计算频率点(i,j)与频域中心的距离
D(i,j)=sqrt((i-round(N1/2)^2+(j-round(N2/2)^2))
for j=1:N2
distance=sqrt((i-n1)^2+(j-n2)^2);
if distance==0
h=0;
else
h=1/(1+(distance/d0)^(2*n));% 根据巴特沃斯低通滤波器公式为
1/(1+[D(i,j)/D0]^2n)
end
s(i,j)=h*s(i,j);% 频域图像乘以滤波器的系数
end
end
% real 函数取元素的实部
s=real(ifft2(ifftshift(s)));% 最后进行二维傅里叶反变换转换为时域图像
subplot(133),imshow(s,[]);
title('巴特沃斯低通滤波图像');