1.频域滤波与图像增强
1.1 基本步骤
(1) zero padding
(2) Spectrum centralization
(3) DFT:
(4) Symmetry
(5)
(6) IDFT :
(7) from with size top-left part.
1.2 MATLAB代码
clc,close all;
I = imread('byh.jpg');% MATLAB库里面有图源 cell.tif
I = rgb2gray(I);
[m,n] = size(I); % m为行，n为列
J = zeros(2*m,2*n);
% 进行扩展
for i=1:m
for j=1:n
J(i,j) = I(i,j);
end
end
% 中心化处理 Spectrum ccentralization
K = J;
for i=1:m
for j=1:n
K(i,j) = I(i,j)*(-1).^(i+j);
end
end
% 傅里叶变换
F = log(1+fft2(im2double(K)));
% H 函数
H = zeros(2*m,2*n);
square = 70; % 滤波函数 H 的边长
H(m-square:m+square,n-square:n+square)=1;
% 相乘
FF = F .* H;
% 反傅里叶变换
ff = ifft2(exp(FF)-1);
ff_r = real(ff);
ff_final = ff_r;
for i=1:2*m
for j=1:2*n
ff_final(i,j) = ff_final(i,j) * (-1)^(i+j);
end
end
subplot(241),imshow(I,[]),xlabel('(a) 原图');
subplot(242),imshow(J,[]),xlabel('(b) 进行扩展操作后的图片');
subplot(243),imshow(K,[]),xlabel('(c) 频域中心化后的图像');
subplot(244),imshow(F,[]),xlabel('(d) 傅里叶变换后的图像');
subplot(245),imshow(H,[]),xlabel('(e) 自制的低通滤波器');
subplot(246),imshow(FF,[]),xlabel('(f) 低通滤波处理后的图像');
subplot(247),imshow(ff_final,[]),xlabel('(g) 反傅里叶变换后的图像');
subplot(248),imshow(ff_final(1:m,1:n),[]),xlabel('(h) 取(g)图中的左上部分');
result:
“如果不进行扩展，直接进行(低通)滤波，结果如下：
clc,close all;
I = imread('byh.jpg');% MATLAB库里面有图源 cell.tif
I = rgb2gray(I);
[m,n] = size(I); % m为行，n为列
J = zeros(m,n);
% 进行扩展
for i=1:m
for j=1:n
J(i,j) = I(i,j);
end
end
% 中心化处理 Spectrum ccentralization
K = J;
for i=1:m
for j=1:n
K(i,j) = I(i,j)*(-1).^(i+j);
end
end
% 傅里叶变换
F = log(1+fft2(im2double(K)));
% H 函数
H = zeros(m,n);
square = 35; % 滤波函数 H 的边长
H(round(m/2)-square:round(m/2)+square,round(n/2)-square:round(n/2)+square)=1;
% 相乘
FF = F .* H;
% 反傅里叶变换
ff = ifft2(exp(FF)-1);
ff_r = real(ff);
ff_final = ff_r;
for i=1:m
for j=1:n
ff_final(i,j) = ff_final(i,j) * (-1)^(i+j);
end
end
subplot(241),imshow(I,[]),xlabel('(a) 原图');
subplot(242),imshow(J,[]),xlabel('(b) 进行扩展操作后的图片');
subplot(243),imshow(K,[]),xlabel('(c) 频域中心化后的图像');
subplot(244),imshow(F,[]),xlabel('(d) 傅里叶变换后的图像');
subplot(245),imshow(H,[]),xlabel('(e) 自制的低通滤波器');
subplot(246),imshow(FF,[]),xlabel('(f) 低通滤波处理后的图像');
subplot(247),imshow(ff_final,[]),xlabel('(g) 反傅里叶变换后的图像');
subplot(248),imshow(ff_final(1:m,1:n),[]),xlabel('(h) 取(g)图中的左上部分');
“只修改了部分代码，表达清楚意思即可hhh ;-)
result:
2 频域滤波器分类
2.1 低通滤波器
“理想低通滤波器(ILPF)特点：物理上不可实现；有抖动现象；滤除高频成分使图像变模糊
其中，
D0为截止频率；P = 2M，Q = 2N，M,N为待滤波图像的尺寸。
clc,close all;
I = imread('kobi.png');% MATLAB库里面有图源 cell.tif byh.jpg
I = rgb2gray(I);
[P,Q] = size(I); % m为行，n为列
F = log(1+fftshift(fft2(I))); % 显示频谱
H = ILPF_ZH(P,Q,70); % 滤波矩阵 / 系统传递函数
F_new = F .* H; % 滤波后的频谱
F_ishift = ifft2(exp(F_new)-1);
ff_r = real(F_ishift);
ff_final = ff_r;
for i=1:P
for j=1:Q
ff_final(i,j) = ff_final(i,j) * (-1)^(i+j);
end
end
imshow(ff_final,[])
subplot(231),imshow(I,[]),xlabel('(a) 原图');
subplot(232),imshow(F,[]),xlabel('(b) 原图的频谱');
subplot(233),imshow(H,[]),xlabel('(c) 滤波器');
subplot(234),x=1:1:P;y=1:1:Q;[X,Y]=meshgrid(x,y);mesh(X,Y,H');xlabel('(d) 滤波器三维展示');
subplot(235),imshow(F_new,[]),xlabel('(e) 滤波后的频谱');
subplot(236),imshow(ff_final,[]),xlabel('(f) 滤波后的图像')
%% 理想低通滤波器
function H = ILPF_ZH(P,Q,D0)
% P = 2551;Q = 2551;% P和Q均为奇数
H = zeros(P,Q);
% D0 = 20;
for i=1:P
for j=1:Q
D = ((i - round((P+1)/2))^2 + (j - round((Q+1)/2))^2)^0.5;
if D H(i,j) = 1;
end
end
end
end
RESULT
“巴特沃思(Butterworth)低通滤波器
BLPF
其中，
clc,close all;
I = imread('kobi.png');% MATLAB库里面有图源 cell.tif byh.jpg
I = rgb2gray(I);
[P,Q] = size(I); % m为行，n为列
F = log(1+fftshift(fft2(I))); % 显示频谱
H = BLPF_ZH(P,Q,70,2); % 巴特沃思(Butterworth)低通滤波矩阵 / 系统传递函数
F_new = F .* H; % 滤波后的频谱
F_ishift = ifft2(exp(F_new)-1);
ff_r = real(F_ishift);
ff_final = ff_r;
for i=1:P
for j=1:Q
ff_final(i,j) = ff_final(i,j) * (-1)^(i+j);
end
end
imshow(ff_final,[])
subplot(231),imshow(I,[]),xlabel('(a) 原图');
subplot(232),imshow(F,[]),xlabel('(b) 原图的频谱');
subplot(233),imshow(H,[]),xlabel('(c) 滤波器');
subplot(234),x=1:1:P;y=1:1:Q;[X,Y]=meshgrid(x,y);mesh(X,Y,H');xlabel('(d) 滤波器三维展示');
subplot(235),imshow(F_new,[]),xlabel('(e) 滤波后的频谱');
subplot(236),imshow(ff_final,[]),xlabel('(f) 滤波后的图像')
%% 巴特沃思(Butterworth)低通滤波器
function H = BLPF_ZH(P,Q,D0,n)
H = zeros(P,Q);
for i=1:P
for j=1:Q
D = ((i - round((P+1)/2))^2 + (j - round((Q+1)/2))^2)^0.5;
H(i,j) = 1/(1+(D/D0)^(2*n));
end
end
end
RESULT
“高斯低通滤波器(GLPF)
其中，
clc,close all;
I = imread('kobi.png');% MATLAB库里面有图源 cell.tif byh.jpg
I = rgb2gray(I);
[P,Q] = size(I); % m为行，n为列
F = log(1+fftshift(fft2(I))); % 显示频谱
H = GLPF_ZH(P,Q,50); % 高斯低通滤波器(GLPF)
F_new = F .* H; % 滤波后的频谱
F_ishift = ifft2(exp(F_new)-1);
ff_r = real(F_ishift);
ff_final = ff_r;
for i=1:P
for j=1:Q
ff_final(i,j) = ff_final(i,j) * (-1)^(i+j);
end
end
imshow(ff_final,[])
subplot(231),imshow(I,[]),xlabel('(a) 原图');
subplot(232),imshow(F,[]),xlabel('(b) 原图的频谱');
subplot(233),imshow(H,[]),xlabel('(c) 滤波器');
subplot(234),x=1:1:P;y=1:1:Q;[X,Y]=meshgrid(x,y);mesh(X,Y,H');xlabel('(d) 滤波器三维展示');
subplot(235),imshow(F_new,[]),xlabel('(e) 滤波后的频谱');
subplot(236),imshow(ff_final,[]),xlabel('(f) 滤波后的图像')
%% 高斯低通滤波器(GLPF)
function H = GLPF_ZH(P,Q,D0)
H = zeros(P,Q);
for i=1:P
for j=1:Q
D = ((i - round((P+1)/2))^2 + (j - round((Q+1)/2))^2)^0.5;
H(i,j) = exp(-(D)^2 / (2 * D0^2));
end
end
end
RESULT
2.2 高通滤波器
通常，频域高通滤波器可以看成是频域低通滤波器的反操作，即
“理想高通滤波器(IHPF)
clc,close all;
I = imread('kobi.png');% MATLAB库里面有图源 cell.tif byh.jpg
I = rgb2gray(I);
[P,Q] = size(I); % m为行，n为列
F = log(1+fftshift(fft2(I))); % 显示频谱
H = IHPF_ZH(P,Q,15); % 理想高通滤波矩阵 / 系统传递函数
F_new = F .* H; % 滤波后的频谱
F_ishift = ifft2(exp(F_new)-1);
ff_r = real(F_ishift);
ff_final = ff_r;
for i=1:P
for j=1:Q
ff_final(i,j) = ff_final(i,j) * (-1)^(i+j);
end
end
imshow(ff_final,[])
subplot(231),imshow(I,[]),xlabel('(a) 原图');
subplot(232),imshow(F,[]),xlabel('(b) 原图的频谱');
subplot(233),imshow(H,[]),xlabel('(c) 滤波器');
subplot(234),x=1:1:P;y=1:1:Q;[X,Y]=meshgrid(x,y);mesh(X,Y,H');xlabel('(d) 滤波器三维展示');
subplot(235),imshow(F_new,[]),xlabel('(e) 滤波后的频谱');
subplot(236),imshow(ff_final,[]),xlabel('(f) 滤波后的图像')
%% 理想高通滤波器
function H = IHPF_ZH(P,Q,D0)
H = ones(P,Q);
for i=1:P
for j=1:Q
D = ((i - round((P+1)/2))^2 + (j - round((Q+1)/2))^2)^0.5;
if D H(i,j) = 0;
end
end
end
end
RESULT
“巴特沃思(Butterworth)高通滤波器
clc,close all;
I = imread('kobi.png');% MATLAB库里面有图源 cell.tif byh.jpg
I = rgb2gray(I);
[P,Q] = size(I); % m为行，n为列
F = log(1+fftshift(fft2(I))); % 显示频谱
H = BHPF_ZH(P,Q,8,2); % 巴特沃思(Butterworth)高通滤波矩阵 / 系统传递函数
F_new = F .* H; % 滤波后的频谱
F_ishift = ifft2(exp(F_new)-1);
ff_r = real(F_ishift);
ff_final = ff_r;
for i=1:P
for j=1:Q
ff_final(i,j) = ff_final(i,j) * (-1)^(i+j);
end
end
imshow(ff_final,[])
subplot(231),imshow(I,[]),xlabel('(a) 原图');
subplot(232),imshow(F,[]),xlabel('(b) 原图的频谱');
subplot(233),imshow(H,[]),xlabel('(c) 滤波器');
subplot(234),x=1:1:P;y=1:1:Q;[X,Y]=meshgrid(x,y);mesh(X,Y,H');xlabel('(d) 滤波器三维展示');
subplot(235),imshow(F_new,[]),xlabel('(e) 滤波后的频谱');
subplot(236),imshow(ff_final,[]),xlabel('(f) 滤波后的图像')
%% 巴特沃思(Butterworth)高通滤波器
function H = BHPF_ZH(P,Q,D0,n)
H = zeros(P,Q);
for i=1:P
for j=1:Q
D = ((i - round((P+1)/2))^2 + (j - round((Q+1)/2))^2)^0.5;
H(i,j) = 1/(1+(D0/D)^(2*n));
end
end
end
RESULT
“高斯高通滤波器(GLPF)
clc,close all;
I = imread('kobi.png');% MATLAB库里面有图源 cell.tif byh.jpg
I = rgb2gray(I);
[P,Q] = size(I); % m为行，n为列
F = log(1+fftshift(fft2(I))); % 显示频谱
H = GHPF_ZH(P,Q,8); % 高斯高通滤波器(GLPF)
F_new = F .* H; % 滤波后的频谱
F_ishift = ifft2(exp(F_new)-1);
ff_r = real(F_ishift);
ff_final = ff_r;
for i=1:P
for j=1:Q
ff_final(i,j) = ff_final(i,j) * (-1)^(i+j);
end
end
imshow(ff_final,[])
subplot(231),imshow(I,[]),xlabel('(a) 原图');
subplot(232),imshow(F,[]),xlabel('(b) 原图的频谱');
subplot(233),imshow(H,[]),xlabel('(c) 滤波器');
subplot(234),x=1:1:P;y=1:1:Q;[X,Y]=meshgrid(x,y);mesh(X,Y,H');xlabel('(d) 滤波器三维展示');
subplot(235),imshow(F_new,[]),xlabel('(e) 滤波后的频谱');
subplot(236),imshow(ff_final,[]),xlabel('(f) 滤波后的图像')
%% 高斯高通滤波器(GHPF)
function H = GHPF_ZH(P,Q,D0)
H = zeros(P,Q);
for i=1:P
for j=1:Q
D = ((i - round((P+1)/2))^2 + (j - round((Q+1)/2))^2)^0.5;
H(i,j) = 1 - exp(-(D)^2 / (2 * D0^2));
end
end
end
RESULT
汇合在一起的代码如下：
clc,close all;
I = imread('kobi.png');% MATLAB库里面有图源 cell.tif byh.jpg
I = rgb2gray(I);
[P,Q] = size(I); % m为行，n为列
F = log(1+fftshift(fft2(I))); % 显示频谱
% H = ILPF_ZH(P,Q,70); % 理想低通滤波矩阵 / 系统传递函数
% H = BLPF_ZH(P,Q,70,2); % 巴特沃思(Butterworth)低通滤波矩阵 / 系统传递函数
% H = GLPF_ZH(P,Q,50); % 高斯低通滤波器(GLPF)
% H = IHPF_ZH(P,Q,15); % 理想高通滤波矩阵 / 系统传递函数
% H = BHPF_ZH(P,Q,8,2); % 巴特沃思(Butterworth)高通滤波矩阵 / 系统传递函数
H = GHPF_ZH(P,Q,8); % 高斯高通滤波器(GLPF)
F_new = F .* H; % 滤波后的频谱
F_ishift = ifft2(exp(F_new)-1);
ff_r = real(F_ishift);
ff_final = ff_r;
for i=1:P
for j=1:Q
ff_final(i,j) = ff_final(i,j) * (-1)^(i+j);
end
end
imshow(ff_final,[])
subplot(231),imshow(I,[]),xlabel('(a) 原图');
subplot(232),imshow(F,[]),xlabel('(b) 原图的频谱');
subplot(233),imshow(H,[]),xlabel('(c) 滤波器');
subplot(234),x=1:1:P;y=1:1:Q;[X,Y]=meshgrid(x,y);mesh(X,Y,H');xlabel('(d) 滤波器三维展示');
subplot(235),imshow(F_new,[]),xlabel('(e) 滤波后的频谱');
subplot(236),imshow(ff_final,[]),xlabel('(f) 滤波后的图像')
%% 理想低通滤波器
function H = ILPF_ZH(P,Q,D0)
H = zeros(P,Q);
for i=1:P
for j=1:Q
D = ((i - round((P+1)/2))^2 + (j - round((Q+1)/2))^2)^0.5;
if D H(i,j) = 1;
end
end
end
end
%% 巴特沃思(Butterworth)低通滤波器
function H = BLPF_ZH(P,Q,D0,n)
H = zeros(P,Q);
for i=1:P
for j=1:Q
D = ((i - round((P+1)/2))^2 + (j - round((Q+1)/2))^2)^0.5;
H(i,j) = 1/(1+(D/D0)^(2*n));
end
end
end
%% 高斯低通滤波器(GLPF)
function H = GLPF_ZH(P,Q,D0)
H = zeros(P,Q);
for i=1:P
for j=1:Q
D = ((i - round((P+1)/2))^2 + (j - round((Q+1)/2))^2)^0.5;
H(i,j) = exp(-(D)^2 / (2 * D0^2));
end
end
end
%% 理想高通滤波器
function H = IHPF_ZH(P,Q,D0)
H = ones(P,Q);
for i=1:P
for j=1:Q
D = ((i - round((P+1)/2))^2 + (j - round((Q+1)/2))^2)^0.5;
if D H(i,j) = 0;
end
end
end
end
%% 巴特沃思(Butterworth)高通滤波器
function H = BHPF_ZH(P,Q,D0,n)
H = zeros(P,Q);
for i=1:P
for j=1:Q
D = ((i - round((P+1)/2))^2 + (j - round((Q+1)/2))^2)^0.5;
H(i,j) = 1/(1+(D0/D)^(2*n));
end
end
end
%% 高斯高通滤波器(GHPF)
function H = GHPF_ZH(P,Q,D0)
H = zeros(P,Q);
for i=1:P
for j=1:Q
D = ((i - round((P+1)/2))^2 + (j - round((Q+1)/2))^2)^0.5;
H(i,j) = 1 - exp(-(D)^2 / (2 * D0^2));
end
end
end
END
