实验二 频率域处理
一.实验目的及要求
1.掌握MATLAB傅里叶变换方法 2.掌握频率域滤波的概念及方法
3.熟练掌握频率域各类滤波器的特点及设计方法 4.利用MATLAB程序进行频域滤波
二、实验内容
㈠ 研究以下程序,分析程序功能;输入执行各命令行,认真观察命令执行的结果。熟悉程序中所使用函数的调用方法,改变有关参数,观察试验结果。 1.简单二值图像的傅里叶变换及幅度谱
clear all; close all; clc; %创建一个二值图像 f = zeros(256,256); f(100:156,118:138) = 1; %计算其傅里叶变换 F = fft2(f);
F1 = fftshift(F); %中心化 %显示原图像及其幅度谱 figure, imshow(f)
figure, imshow(abs(F),[]); colormap(gray); colorbar; %将幅度谱对数变换后再显示
figure, imshow(log(1+ abs(F)),[]); colormap(gray); colorbar; figure, imshow(log(1+ abs(F1)),[]); colormap(gray); colorbar; figure, mesh(1:256,1:256, abs(F1)); %将图像旋转45度,然后重做上述步骤 f45 =imrotate(f,45,'crop'); F2 = fft2(f45);
9
F3 = fftshift(F2); figure, imshow(f45);
figure, imshow(log(1+ abs(F3)),[]); colormap(gray); colorbar figure, mesh(1:256,1:256, abs(F3));
%------------------------------------------------------------
2. 频率域平滑(低通)滤波
⑴ 空间域高斯低通滤波器的频率域实现
close all; clear all; clc; I = imread('text.png'); f = im2double(I);
%定义高斯低通滤波器标准差变量,初始化为1 sigma = 2;
%根据标准差sigma大小,确定滤波模板的尺寸 hsize = floor(6*sigma); if mod(hsize,2) == 0 hsize = hsize+1; end
%生成指定大小的空间域高斯低通滤波模板 h = fspecial('gaussian', hsize, sigma); AB = size(f); CD = size(h); PQ = AB + CD -1; F = fft2(f,PQ(1),PQ(2)); H = fft2 (h,PQ(1),PQ(2)); g = real(ifft2(F .* H));
g = imcrop(g, [1,1,AB(2)-1, AB(1)-1]); g =im2uint8(g);
figure, imshow(I); title('Original Image'); figure, imshow(g); title('the Smoothed Image'); %------------------------------------------------------------
⑵ 理想低通滤波器的简单实现
close all; clear all; clc;
I = imread('Fig0333(a)(test_pattern_blurring_orig).tif'); f = im2double(I);
10
D0 = 10; % 定义滤波器截止频率 F = fft2(f); %快速傅里叶变换 Fc = fftshift(F); % 中心化 [height, width]=size(Fc);
%初始化理想低通滤波器系数矩阵(频率域设计滤波器) H = zeros(height, width); %计算理想低通滤波器系数 for u=1:height for v =1:width
if(sqrt((u-height/2)^2 + (v-width/2)^2) < D0) H(u,v) =1; else
H(u,v) = 0;
end end end
%进行理想低通滤波,图像傅里叶变换矩阵与滤波器系数矩阵对应点乘G = Fc.* H
G = ifftshift(G); % 去中心化
g = real(ifft2(G)); % 逆傅里叶变换并取实部
g = im2uint8(g); % 把输出图像转换为uint8数据类型 figure, imshow(I);title('Original Image');
figure, imshow(g); title('the ILP Smoothed Image'); %------------------------------------------------------------
⑶ 巴特沃斯Butterworth低通滤波器的简单实现
close all; clear all; clc;
I = imread('Fig0333(a)(test_pattern_blurring_orig).tif'); f = im2double(I);
D0 = 10; %定义滤波器截止频率,可令D0分别等于 30,60..观察结果 n = 2; % 定义滤波器阶数,可增大阶数n=3,4…观察滤波效果 F = fft2(f); Fc = fftshift(F);
[height, width] = size(Fc);
%初始化巴特沃斯低通滤波器系数矩阵(频率域设计滤波器) H = zeros(height, width);
11
%计算巴特沃斯低通滤波器系数矩阵 for u =1:height
for v =1:width
%计算巴特沃斯低通滤波器在频率点 (u,v)的系数值 Dist = sqrt((u-height/2)^2 + (v-width/2)^2); H(u, v) = 1/(1+(Dist/D0).^(2*n)); end end
%进行巴特沃斯低通滤波, 图像傅里叶变换矩阵与滤波器系数矩阵对应点乘G = Fc.* H; G = ifftshift(G); g = real(ifft2(G)); g = im2uint8(g);
figure, imshow(I);title('Original Image');
figure, imshow(g); title('the Butterworth LPF Smoothed Image'); %------------------------------------------------------------
⑷ 高斯Gaussian低通滤波器的简单实现
close all; clear all; clc;
I = imread('Fig0333(a)(test_pattern_blurring_orig).tif'); f = im2double(I);
D0 = 10; %定义滤波器截止频率,可令D0分别等于 30,60..观察结果 F = fft2(f); Fc = fftshift(F); [height,width]=size(Fc);
%初始化高斯低通滤波器系数矩阵(频率域设计滤波器) H = zeros(height, width); %计算高斯低通滤波器系数矩阵 for u =1:height for v =1:width
%计算高斯高通滤波器在频率点 (u,v)的系数值 sqDist = (u-height/2)^2 + (v-width/2)^2; H(u, v) = exp(-sqDist/(2*D0^2)); end end
%进行高斯低通滤波, 图像傅里叶变换矩阵与滤波器系数矩阵对应点乘 G = Fc.* H;
12
G = ifftshift(G); g = real(ifft2(G)); g = im2uint8(g);
figure, imshow(I);title('Original Image');
figure, imshow(g); title('the Gaussian HPF Smoothed Image'); %------------------------------------------------------------
3. 频率域锐化(高通)滤波
⑴ 高斯Gaussian高通滤波器的简单实现
close all; clear all; clc;
I = imread('Fig0333(a)(test_pattern_blurring_orig).tif'); f = im2double(I);
D0 = 10; %定义滤波器截止频率,可令D0分别等于 30,60..观察结果 F = fft2(f); Fc = fftshift(F); [height,width]=size(Fc);
%初始化高斯高通滤波器系数矩阵(频率域设计滤波器) H = zeros(height, width); %计算高斯高通滤波器系数矩阵 for u =1:height for v =1:width
%计算高斯高通滤波器在频率点 (u,v)的系数值(注意与高斯低通滤波器之间的差异) sqDist = (u-height/2)^2 + (v-width/2)^2; H(u, v) = 1- exp(-sqDist/(2*D0^2)); end end
%进行高斯高低通滤波, 图像傅里叶变换矩阵与滤波器系数矩阵对应点乘 G = Fc.* H; G = ifftshift(G); g = real(ifft2(G)); g = im2uint8(g);
figure, imshow(I);title('Original Image');
figure, imshow(g); title('the Gaussian HPF Image'); %------------------------------------------------------------
⑵ Laplacian 图像锐化增强的频域实现
clear all; close all; clc;
13