前言：作为一名图像算法攻城狮，那是在2014年大三下学期，一本深绿色的《数字图像处理》（冈萨雷斯的英文版）出现在自己课桌前。偶然间打开的一扇意外之门，就这样结下了一段不解之缘，那些日子不断上网搜代码的自己，那个刚上机不到二十分钟就把作业提交的自己，早已随往日的岁月飘扬而去。三年的烟酒僧，两年的酱油工，而今只觉脑子越来越不够用，这次决心回炉重造，希望能够通过固本培基，打开思路，话不多说，开始上课！ ----2020-7-6
一、实验目的 通过本实验使学生掌握使用MATLAB的二维傅里叶变换进行频域增强的方法。 二、实验原理 本实验是基于数字图像处理课程中的图像频域增强理论来设计的。 本实验的准备知识：第四章 频域图像增强中的一维傅里叶变换和二维傅里叶变换，频域图像增强的步骤，频域滤波器。根据教材285页到320页的内容，开展本实验。 可能用到的函数： 1、 延拓函数 padarray 例：A=[1,2;3,4]; B=padarray(A,[2,3],’post’); 则结果为 B = 1 2 0 0 0 3 4 0 0 0 0 0 0 0 0 0 0 0 0 0 使用该函数实现图像的0延拓。Padarray还有其它用法，请用help查询。
2 、低通滤波器生成函数 首先编写dftuv函数，如下 function [U,V]=dftuv(M,N) %DFTUV Computes meshgrid frequency matrices. % [U,V]=DFTUV(M,N] computes meshgrid frequency matrices U and V. U and V are useful for computing frequency-domain filter functions that can be used with DFTFILT. U and V are both M-by-N.
% Set up range of variables. u=0:(M-1); v=0:(N-1); % Compute the indices for use in meshgrid. idx=find(u>M/2); u(idx)=u(idx)-M; idy=find(v>N/2); v(idy)=v(idy)-N; %Compute the meshgrid arrays. [V,U]=meshgrid(v,u); 然后编写低通滤波器函数 function [H,D]=lpfilter(type,M,N,D0,n) % LPFILTER computers frequency domain lowpass filters. % H=lpfilter(TYPE,M,N,D0,n) creates the transfer function of a lowpass filter, H, of the specified TYPE and size(M-by-N). To view the filter as an image or mesh plot, it should be centered using H=fftshift(H). % valid values for TYPE, D0, and n are: % ‘ideal’ Ideal lowpass filter with cutoff frequency D0. n need not be supplied. D0 must be positive. % ‘btw’ Butterworth lowpass filter of ordern, and cutoff D0. The default value for n is 1. D0 must be positive. % ‘gaussian’ Gaussian lowpass filter with cutoff (standard deviation)D0. n need not be supplied. D0 must be positive. %Use function dftuv to set up the meshgrid arrays needed for computing the required distances. [U,V]=dftuv(M,N); % D=sqrt(U.2+V.2); % Compute the distances D(U,V) % Begin filter computations. switch type case ‘ideal’ H=double(D<=D0); case ‘btw’ if nargin==4 n=1; end H=1./(1+(D./D0).^(2*n)); case ‘gaussian’ H=exp(-(D.2)./(2*(D02))); otherwise error(‘Unknown filter type’) end 通过调用函数lpfilter可生成相应的滤波器掩膜矩阵。 参考该函数可相应的生成高通滤波器函数。
3、 频域滤波 F=fft2(f,size(H,1),size(H,2)); % 对延拓的 f 计算 FFT。注意，这里隐含着对 f 的延拓。 G=real(ifft2(H.*F)); % 滤波 Gf=G(1:size(f,1),1:size(f,2)); %裁剪后的图像
三、实验内容 （一）图像频域增强的步骤 参考教材286页的Figure 4.36，重复该图像中的步骤，并将相应的结果显示出来。
（二）频域低通滤波 产生实验四中的白条图像。设计不同截止频率的理想低通滤波器、Butterworth低通滤 波器，对其进行频域增强。观察频域滤波效果，并解释之。
程序1：
I=zeros(64,64); I(32-20:32+20,32-8:32+8)=255; subplot(1,2,1) imshow(I); s=fftshift(fft2(I)); [M,N]=size(s); n1=floor(M/2); n2=floor(N/2); d0=5; for i=1:M for j=1:N d=sqrt((i-n1)^2+(j-n2)^2); h=1*exp(-1/2*(d^2/d0^2)); s(i,j)=h*s(i,j); end end s=ifftshift(s); s=uint8(real(ifft2(s))); subplot(1,2,2); imshow(s)程序2：
I=zeros(64,64); I(32-20:32+20,32-8:32+8)=255; subplot(1,2,1) imshow(I); s=fftshift(fft2(I)); [M,N]=size(s); n1=floor(M/2); n2=floor(N/2); d0=50; for i=1:M for j=1:N d=sqrt((i-n1)^2+(j-n2)^2); h=1*exp(-1/2*(d^2/d0^2)); s(i,j)=h*s(i,j); end end s=ifftshift(s); s=uint8(real(ifft2(s))); subplot(1,2,2); imshow(s);程序3：
I=zeros(64,64); I(32-20:32+20,32-8:32+8)=255; subplot(1,2,1) imshow(I); s=fftshift(fft2(I)); [M,N]=size(s); n1=floor(M/2); n2=floor(N/2); d0=5; for i=1:M for j=1:N d=sqrt((i-n1)^2+(j-n2)^2); h=1/(1+0.414*(d/d0)^(2*n1)); s(i,j)=h*s(i,j); end end s=ifftshift(s); s=uint8(real(ifft2(s))); subplot(1,2,2); imshow(s);设计不同截止频率的理想低通滤波器、Butterworth低通滤波器，对含高斯噪声的lena图像进行频域增强。观察频域滤波效果，并解释之。
程序1：
clear all; A=imread('F:\Images\lena.bmp'); subplot(1,3,1),imshow(A,[]); title('原图') B=imnoise(A,'gauss',0.02); subplot(1,3,2),imshow(B,[]); title('添加高斯噪声后的图像') f=double(B); g=fft2(f); g=fftshift(g); [N1,N2]=size(g); n=2; d0=50; n1=fix(N1/2); n2=fix(N2/2); for i=1:N1 for j=1:N2 d=sqrt((i-n1)^2+(j-n2)^2); c=double(d<=d0); result(i,j)=c*g(i,j); end end resul=ifftshift(result); X2=ifft2(result); X3=uint8(real(X2)); subplot(1,3,3),imshow(X3); title('D0=50,低通滤波器')程序2： Butterworth低通滤波器程序：
clc,clear,close all; I=imread('F:\Images\lena.bmp'); %原图为灰度图像 II=imnoise(I,'gauss',0.02); I1=medfilt2(II,[3,3]); subplot(1,3,1) imshow(I,[]);%把图像显示出来 subplot(1,3,2) imshow(I1,[]); f=double(I1); %图像存储类型转换 g=fft2(f); %傅立叶变换 g=fftshift(g); %转换数据矩阵 [N1,N2]=size(g); %测量图像尺寸参数 n=2; d0=50; n1=fix(N1/2); n2=fix(N2/2); for i=1:N1 for j=1:N2 d=sqrt((i-n1)^2+(j-n2)^2); h=1/(1+0.414*(d/d0)^(2*n));%计算Butterworth低通转换函数 result(i,j)=h*g(i,j); end end result=ifftshift(result); X2=ifft2(result); X3=uint8(real(X2)); subplot(1,3,3) imshow(X3)%显示频域增强后的图像（三） 频域高通滤波 设计不同截止频率的理想高通滤波器、Butterworth高通滤波器，对上述白条图像进行频域增强。观察频域滤波效果，并解释之。 程序1：
I=zeros(64,64); I(32-20:32+20,32-8:32+8)=255; subplot(1,2,1); imshow(I); s=fftshift(fft2(I)); [M,N]=size(s); n1=floor(M/2); n2=floor(N/2); d0=5; for i=1:M for j=1:N d=sqrt((i-n1)^2+(j-n2)^2); h=1-(1*exp(-1/2*(d^2/d0^2))); s(i,j)=h*s(i,j); end end s=ifftshift(s); s=uint8(real(ifft2(s))); subplot(1,2,2); imshow(s);程序2：
I=zeros(64,64); I(32-20:32+20,32-8:32+8)=255; subplot(1,2,1); imshow(I); s=fftshift(fft2(I)); [M,N]=size(s); n1=floor(M/2); n2=floor(N/2); d0=5; for i=1:M for j=1:N d=sqrt((i-n1)^2+(j-n2)^2); h=1-(1/(1+0.414*(d/d0)^(2*n))); s(i,j)=h*s(i,j); end end s=ifftshift(s); s=uint8(real(ifft2(s))); subplot(1,2,2); imshow(s)设计不同截止频率的理想高通滤波器、Butterworth高通滤波器，对含高斯噪声的lena图像进行频域增强。观察频域滤波效果，并解释之。
程序1：
I1=imread('F:\Images\lena.bmp'); I=imnoise(I1, 'gauss', 0.02); subplot(1,2,1); imshow(I); s=fftshift(fft2(I)); [M,N]=size(s); n1=floor(M/2); n2=floor(N/2); d0=5; for i=1:M for j=1:N d=sqrt((i-n1)^2+(j-n2)^2); h=1-(1*exp(-1/2*(d^2/d0^2))); s(i,j)=h*s(i,j); end end s=ifftshift(s); s=uint8(real(ifft2(s))); subplot(1,2,2); imshow(s);程序2：
I1=imread('F:\Images\lena.bmp'); I=imnoise(I1, 'gauss', 0.02); subplot(1,2,1); imshow(I); s=fftshift(fft2(I)); [M,N]=size(s); n1=floor(M/2); n2=floor(N/2); d0=5; for i=1:M for j=1:N d=sqrt((i-n1)^2+(j-n2)^2); h=1-(1/(1+0.414*(d/d0)^(2*n))); s(i,j)=h*s(i,j); end end s=ifftshift(s); s=uint8(real(ifft2(s))); subplot(1,2,2); imshow(s);**四、实验方法与步骤 1、顺序完成上述实验内容 2、按照实验内容要求，分析编程，将程序和实验结果整理成word文档，分析结果，编写实验报告。
五、实验报告要求 1、本实验由学生单人独立完成。 2、每个实验均按统一格式编写实验报告。 实验报告内容包括：实验要求，实验项目，典型程序流程图，程序清单，数据结果和分析讨论。**
