数字图像处理指导-测控10版本.doc

数字图像处理指导-测控10版本 作者： 日期：37 个人收集整理 勿做商业用途 实验一：图像预处理实验（验证性实验）实验目的：要求学生掌握图像读取、显示和保存的方法，掌握空间域增强方法和频率域增强方法.实验要求：上机运行,调试通过。实验内容：1） 对附录一中的练习111、练习13，进行上机运行和调试.2） 对上述练习，改变一些函数的参数，观察运行结果的变化。提交试验报告的内容：1、 实现对一幅图（比如：pout.tif）进行直方图均衡化的源代码；2、 显示直方图均衡化前后的图像及其直方图；3、 针对处理前后的图像和对应的直方图进行分析,得出直方图均衡化的作用。实验二、图像的傅立叶变换（验证性实验)实验目的1了解图像变换的意义和手段;2熟悉傅里叶变换的基本性质;3熟练掌握FFT的方法与应用；4理解频率域高低频与时域图像的对应关系;5通过本实验掌握利用MATLAB编程实现数字图像的傅立叶变换。实验原理1、应用傅立叶变换进行图像处理傅里叶变换是线性系统分析的一个有力工具，它能够定量地分析诸如数字化系统、采样点、电子放大器、卷积滤波器、噪音和显示点等的作用。通过实验培养这项技能，将有助于解决大多数图像处理问题。对任何想在工作中有效应用数字图像处理技术的人来说,把时间用在学习和掌握博里叶变换上是很有必要的。2、傅立叶（Fourier）变换的定义对于二维信号，二维Fourier变换定义为：二维离散傅立叶变换为: 图像的傅立叶变换与一维信号的傅立叶变换一样,有快速算法,具体参见参考书目，有关傅立叶变换的快速算法的程序不难找到.实际上,现在有实现傅立叶变换的芯片，可以实时实现傅立叶变换.3、用matlab编程观察图像频域中高低频率与时域图像的对应关系：-% 本代码主要演示图像fft2变换的频谱图中高频和低频分别对应 图像时域图中的哪些部分,加深对频域的理解-clear；close all;A=imread(lena.jpg)； %读入图像A=rgb2gray(A);figure（1);imshow（A）; %显示图像A=double（A）;%- 对图像进行fft变换，看看频谱图和ifft变换的结果 -F0=fft2（A)； 计算二维傅里叶变换F1=fftshift(F0); %将直流分量移到频谱图的中心figure（2)， imshow(log（abs(F1)+1),0 10）；显示变换后的频谱图F0=fftshift（F1)；R0=ifft2(F0）； %进行反变换figure(3）；imshow（uint8（R0）); %将变换的结果转换为uint8数据类型,并显示，比较反变换的结果与原图。%- 把频谱中直流和低频系数变为0，只保留高频系数，然后反变化，% 看看高频系数对应图像时域中的哪一部分？%-m，n=size（R0）; 得到图像大小i0=fix(m/2); %图像中心点的坐标j0=fix(n/2）；F2=F1;r=50; 设置系数改变区域的半径。%把中心点附近2W2W的低频全部设置为0；for i=1:m for j=1：n d=sqrt(i-i0)2+(j-j0)2）； if（d=r） F2（i,j)=0; end; end;end；F3=fftshift（F2）; 将直流分量移到频谱图的四个角R1=ifft2（F3）; %进行反变换figure（4）， imshow（log（abs(F2）+1）,0 10）;显示去除低频的频谱图figure(5）；imshow（uint8(R1）); 将变换的结果转换为uint8数据类型,并显示，比较反变换的结果与原图。%- 把频谱中高频系数变为0,只保留低频系数，然后反变化，% 看看低频系数对应图像时域中的哪一部分？-m,n=size(R0）; 得到图像大小i0=fix（m/2）； %图像中心点的坐标j0=fix(n/2);F5=F1;r=50; %设置系数不改变区域的半径。把中心点附近2W*2W的低频全部设置为0；for i=1：m for j=1：n d=sqrt（i-i0）2+（j-j0）2）； if(dr） F5（i,j）=0； end; end;end;F6=fftshift(F5); 将直流分量移到频谱图的四个角R2=ifft2(F6)； %进行反变换figure(6), imshow（log（abs(F5）+1）,0 10)；显示去除低频的频谱图figure（7);imshow(uint8(R2）)； 将变换的结果转换为uint8数据类型，并显示，比较反变换的结果与原图。%- 把高频和低频的反变换结果相加，与原图相比较，看看结果如何？%-R3=R1+R2;figure（8）；imshow(uint8（R3）)；实验内容1打开计算机，安装和启动MATLAB程序；程序组中“work”文件夹中应有待处理的图像文件;2分析、运行其中的代码,观察结果。3 a)。对三幅图像内容区别较大图做fft并显示其频谱； b) 讨论时域图像内容与FFT频谱之间的对应关系。4记录和整理实验报告。提交试验报告的内容：1、对图像内容区别较大的三幅图分别做fft并显示其频谱.2、讨论不同的图像内容的FFT频谱图之间的区别。实验三、数字水印（设计性实验）实验要求1了解数字水印的基本技术；2设计一个基于LSB的数字水印嵌入，提取、攻击方案；3编写代码实现基于LSB的数字水印嵌入、提取功能；参考资料：1、 matlab编程资料.2、 黄仿元，基于LSB的数字水印算法及MATLAB实现。提交的报告内容：1、 数字水印的基本知识介绍。2、 水印嵌入、提取设计方案3、 LSB数字水印的嵌入、提取代码及结果图。4、水印图像被攻击后提取的水印结果图。 附录一：MATLAB数字图像处理1 MATLAB简介11 主要用途及特点主要用途：算法研究主要特点：语句功能强大;不能生成可执行文件.12 MATLAB工作环境 121 Matlab桌面 桌面包括5个子窗口：命令窗口、工作空间窗口、当前目录窗口、历史命令窗口、一个或多个图形窗口（仅在用户显示图形式出现）。命令窗口是用户在提示符（)处键入MATLAB命令和表达式的地方,也是显示那些命令输出的地方。工作空间窗口显示对话中创建的变量和它们的某些信息。当前目录窗口显示当前目录的内容（即路径）。历史命令窗口包含用户已在命令窗口中输入的命令的纪录. 122 使用MATLAB编辑器创建M文件 123 获得帮助（1）help 函数名（2）look 关键词2 数字图像的读取、显示、保存、数据类型和图像类型、数据类型间的转换、图像类型间的转换21 图像的读取语法：imread ( filename ) 说明：读取图像格式名称 描述 可识别扩展符TIFF加标记的图像文件格式.tif， 。tiffJPEG联合图像专家组。jpg, 。jpegGIF图形交换格式.gifBMPWindows位图.bmpPNG可移植网络图形.pngXWDX Window 转储。xwd 语法:M, N=size （ filename ) 说明：给出一幅图像的行数和列数22 图像的显示语法：imshow （ f， G） imshow （f, low high） imshow (f, )说明：G是显示该图像的灰度级数； 小于或等于low的值都显示为黑色，大于或等于high的值都显示为白色. 自动将变量low设置为f的最小值,将high设置为f的最大值.23 图像的保存语法：imwrite （ f, filename） 说明：练习1 f = imread ( saturn。tif)； M, N = size( f ）; g = imread （ trees.tif ）； imshow （ f ); figure， imshow (g); 显示另一幅图像 imwrite （f, s2.jpg)24 数据类型名称 描述double双精度浮点数，范围为10exp（308) 10exp(308)，8字节uint8无符号8比特整数，1字节uint16无符号16比特整数，2字节uint32无符号32比特整数，4字节int8有符号8比特整数，1字节int16有符号16比特整数，2字节int32有符号32比特整数，4字节single单精度浮点数,范围为10exp(38) 10exp(38） ,4字节char字符logical值为0或1四种常用类型：double，uint8,char,logical。25 图像类型名称 描述亮度图像uint8类范围为0 255、double类归一化为0 1二值图像图像取值只有0和1的逻辑数组索引图像索引图像RGB图像彩色图像2 6 数据类与图像类型间的转换261 数据类间的转换语法：B = data_class_name （ A ）举例：若A是8位图像，则B = double ( A ）转换为双精度图像.262 图像类型间的转换函数名称将输入转换为有效的输入图像数据类im2uint8uint8Logical， uint8, uint16和doubleim2uint16uint16Logical, uint8, uint16和doublemat2graydoubledoubleim2doubledoubleLogical， uint8， uint16和doubleim2bwlogicaluint8， uint16和double 练习2f=1 2; 3 4；g=mat2gray(f);gb=im2bw（g, 0.6) %0。6 is a threshold3 亮度变换与空间滤波31 亮度变换函数311 基本亮度变换函数 语法：g=imadjust（f， low_in high_in， low_out high_out，gamma)说明：将图像f中的亮度值影响到g中的新值，即将low_in至high_in之间的值映射到low_out至high_out之间的值， low_in以下的值映射为low_out, high_in以上的值映射为high_out.参数gamma指定了映射曲线的形状。练习3f = imread ( pout。tif)； imshow（f）;g1=imadjust（f， 0 1, 1 0); figure， imshow(g1)；g2=imadjust(f， 0。5 0.75, 1 0， 0。5)； figure, imshow（g2)；312 对比度拉伸变换对数变换方法：g=c*log(1+double（f）)说明：对数变换的一项主要应用是压缩动态范围,常用于傅里叶频谱显示。练习4f = imread （ pout.tif)； F=fft2(f)； Fourier TransformFC=fftshift(F）； %将变换原点移到频率矩形的中心。imshow(abs(FC)， )；S2= log（1+ abs（FC)）；figure，imshow(S2， ）;对比度拉伸方法:g=1。/(1+(m。/（double（f)+eps).E)说明：eps是浮点数的相对精度,可避免f出现0值时的溢出现象。32 直方图处理与函数绘图321 生成并绘制图像的直方图 语法:h=imhist(f， b） p= imhist(f, b）/num1（f) 说明：b是用于形成直方图的灰度级的个数.P是归一化直方图.练习5f = imread ( pout。tif)； imhist（f）；322 直方图均衡化语法：g=histeq(f， nlev）说明：nlev是为输出图像制定的灰度级数。练习6f = imread （ tire。tif）; imshow（f）； figure，imhist（f）; ylim(auto）； %自动设定y轴坐标范围和刻度 g=histeq（f,256）; figure， imshow（g）； figure,imhist（g）； ylim（auto)；33 空间滤波331 线性空间滤波 语法：g=imfilter（f, w， filtering_mode, boundary_options, size_options）说明：w为滤波掩模选项描述滤波类型corr滤波器通过使用相关来完成.该值是默认值。conv滤波器通过使用卷积来完成边界选项P输入图像的边界通过用值P来扩展。P的默认值为0。replicate图像大小通过复制外边界的值来扩展。symmetric图像大小通过反射其边界来扩展。circular图像大小通过将图像看成是一个二维周期函数的一个周期来扩展大小选项full输出图像的大小与被扩展图像的大小相同same输出图像的大小与输入图像的大小相同练习7f = imread （ saturn。tif）； w=ones（31）; 单位矩阵掩模gd=imfilter（f，w);imshow(gd，);332 非线性空间滤波略.34 图像处理工具箱的标准空间滤波器341 线性空间滤波器语法:w=fspecial(type,parameter)说明：见下表练习8f = imread （ saturn。tif）； w4=fspecial(laplacian， 0）w8=1 1 1; 1 8 1; 1 1 1;g4=fimfilter(f， w4， replicate);g8=fimfilter(f， w8, replicate);imshow(f)；figure， imshow（g4）; figure， imshow（g8)； 342 非线性空间滤波器语法：g=ordfilt2（f， order, domain)说明：使用邻域的一组排序元素中的第oder个元素来代替f中的每个元素，而该邻域则由domain中的非零元素指定.语法：g=medfilt2(f， m n，padopt）说明：中值滤波器。padopt指定了三个可能的边界填充选项之一。 zeros(默认值）。 symmetric，镜像反射。 indexed,double类图像,以1填充，否则以0填充。练习9f = imread ( eight.tif）； fn=imnoise(f， saltpepper, 0.2）； 加椒盐噪声gm=medfilt2(fn)；imshow（fn)；figure，imshow（gm）；语法：J = wiener2(I，m n，noise)说明:进行二维适应性去噪过滤处理(维纳滤波）.练习10I = imread（saturn.tif）；J = imnoise（I,gaussian，0，0。005）； 加高斯噪声K = wiener2（J,5 5）；imshow（J)figure， imshow（K)4 频域处理41 在MATLAB中计算并可视化二维DFT练习11 f = imread （ saturn.tif); F=fft2（f)； % Fourier TransformS=abs（F）; %计算傅里叶频谱imshow(S,）；FC=fftshift(F）; 将变换原点移到频率矩形的中心。figure，imshow(abs（FC）， ）；S2= log（1+ abs(FC);figure,imshow（S2, )；频域滤波421 注意事项假设函数f(x，y）和h(x，y）的大小分别为AB和CD，为了避免卷积运算造成的折叠误差（即混叠失真），必须通过对f和g补零，构造两个大小均为PQ的扩充寒暑，才能在频率域中进行运算.其中，P、Q必须满足：422 DFT滤波的基本步骤1、 对f和g进行补零。2、 对f和g分别进行傅里叶变换,得到频域F和H.3、 将变换乘以滤波函数G=H.F；4、 获得G的傅里叶逆变换的实部: g=real(ifft2(G）)；5、 讲左上部的矩形修建为原始大小： g=g（1:size（f，1）,1:size（f，2)）； 练习12：根据DFT滤波的基本步骤，进行编程实现。从空间滤波器获得频率滤波器语法：H=fft2(h， M， N）说明:M、N是滤波器的行数和列数,由被滤波的图像大小决定,是补零的结果。 语法：H=freqz2(h， R, C） 说明：计算FIR滤波器的频率响应。练习13h =1 1 1； 1 -8 1; 1 1 1;H=freqz2（h, 50, 50）；在频率域中直接生成滤波器练习题14:编程实现在频率域中生成巴特沃兹低通滤波器（或高通滤波器，或高通强调滤波器)。5 形态学图像处理51 膨胀和腐蚀511 膨胀语法：A2=imdilate(A, B） 说明：B是结构元（由0和1组成的矩阵）.练习15A=imread(broken_text.tif)；B=0 1 0; 1 1 1; 0 1 0;A2=imdilate（A, B）；imshow（A）;figure,imshow（A2）;512 结构元的生成语法：se=strel（shape, parameter)说明:见下表。513 腐蚀语法：A2=imerode（A, B）练习16A=imread(wirebond_mask.tif);se=strel（disk,10）；A2=imerode(A, se);imshow（A）； figure，imshow（A2）;52 膨胀和腐蚀的组合521 开运算和闭运算语法：C=imopen(A, B) 和 C=imclose（A,B)练习17A=imread（shapes。tif)；se=strel(square，20);fo=imopen（f,se）;imshow（A）; figure，imshow(fo)；fc=imclose(f,se）； figure，imshow（fc)；522 击中击不中变换语法：C=bwhitmiss(A, B1， B2)练习18A=imread(wirebond_mask。tif)；B1=strel(0 0 0 0； 0 1 1 1； 0 1 1 0);B2=strel(1 1 1 1; 1 0 0 0 ; 1 00 0)；g=bwhitmiss（f, B1, B2）； imshow(g）；523 形态学函数语法:g=bwmorph(f, operation, n）说明：该函数可基于膨胀、腐蚀和查找表的组合实现许多有用的操作,n是一个用于指定将被重复操作次数的正整数。Operation说明见下表.练习19A=imread（fingerprint_cleaned.tif）;G1=bwmorph(f, thin, 1）；G2=bwmorph（f， thin, 1);Ginf= bwmorph(f, thin， Inf)；Fs= bwmorph(f， tskel， Inf）;imshow（A); figure,imshow(G1）; figure,imshow（G2)；figure，imshow（Ginf）；figure，imshow（Fs)；53 标记连通分量语法：L, num=bwlabel（f, conn)说明:conn用于指定连接类型(4连接获8连接），num给出所找到的连接分量总数，L成为标记矩阵。练习20A=imread（objects。tif);L，n=bwlabel(f）;r,c=find（L= =3）;rbar=mean(r）;cbar=mean（c)；plot(cbar， rbar， Marker, o, MarkerEdgeColor, k, MarkerFaceColor, k，MarkerSizer，10）;plot（cbar, rbar, Marker， ， MarkerEdgeColor， w)；54 形态学重构语法：out=imreconstruct(marker， mask）541 由重构做开运算练习21f=imread(book_text_bw.tif）;fe=imerode(f,ones(51，1）；fo=imopen（f，ones（51,1）);fobr=imreconstruct（fe, f);imshow(f）;figure，imshow（fe)；figure,imshow（fo）；figure,imshow(fobr);541 填充孔洞 语法：g=imfill（f, holes)542 清除边界对象 语法:g=imclearborder(f, conn）6 图像分割61 点、线和边缘检测611 点检测 点检测模板w:-11-1-18-1-1-1-1检测方法： g=abs(imfilter（double(f)， w)）=T练习22f=imread（moon.tif);w=-1 -1 -1; -1 8 -1; 1 -1 1；g=abs(imfilter(double(f）， w）；T=max(g(：）)；T=T0。9;g=g=T；imshow（f); figure， imshow（g）;612 线检测水平模板、+45度模板、垂直模板、45度模板。练习23f=imread（wirebond_mask。tif);imshow(f）;w=2 1 -1; -1 2 -1; -1 1 2;g=abs(imfilter(double(f）， w)；figure,imshow(g）;613 使用edge函数的边缘检测语法：g，t=edge(f， method， parameter）说明：g是一个逻辑数组，其值为：在f中检测到边缘的位置为1，其他位置为零；t是edge是用的阈值；method为边缘监测器方法，可选为： sobel, prewit, roberts, log(LoG）, zerocoss， canny等；parameter包含两部分：T为指定的阈值，第二部分为dir(检测边缘的首选方向： horizontal, vertical， both），或sigma（标准方差），或H（指定的滤波函数）.练习24f=imread（rice。tif)；imshow（f)；gsobel,t=edge(f， sobel）；figure， imshow（gsobel）；glog,t=edge（f， log);figure， imshow(glog）;gcanny,t=edge(f， canny);figure, imshow(gcanny);62 使用Hough变换的线检测练习25设计与实现一个基于Hough变换的直线检测器.63 阈值处理631 全局阈值处理语法:T=graythresh(f）说明:T是阈值，归一化为0至1之间的值.632 局部阈值处理通过一个形态学顶帽算子并对得到的结果使用graythresh来计算.练习26f=imread(cell。tif);imshow(f);T=graythresh(f）;G=f=T；figure， imshow(g)；64 基于区域的分割641 区域生长 642 区域分裂和合并练习27设计与实现一个基于区域生长的分割程序.65 使用分水岭变换的分割练习28f=imread（cell。tif);imshow(f）;g=im2bw(f, graythresh（f）;figure,imshow(g)；gc=g;D=bwdist（gc）;L=watershed(-D);w=L= =0；g2=g&w；figue，imshow(g2）；66 分割后处理语法：BW2 = bwfill(BW1,c，r,n)说明：填充二进制图像的背景色。（形态学处理）练习29BW1 =1 0 0 0 0 0 0 01 1 1 1 1 0 0 01 0 0 0 1 0 1 01 0 0 0 1 1 1 01 1 1 1 0 1 1 11 0 0 1 1 0 1 01 0 0 0 1 0 1 01 0 0 0 1 1 1 0BW2 = bwfill（BW1,3，3,8）I = imread(blood1.tif)；BW3 = im2bw(I）；BW4 = bwfill(BW3，holes);imshow(BW3）figure， imshow(BW4）语法：bwareaopen说明:二进制图像区域打开，清除小物体。7 表示与描述71 二进制图像对象的参数测量语法：total = bwarea(BW）说明：计算二进制图像对象的面积。语法:eul = bweuler(BW,n）说明：计算二进制图像的欧拉数.练习30BW = imread(circles.tif）;imshow(BW);total=bwarea（BW）； total =15799eul=bweuler（BW）; % eul =2 语法:BW2 = bwperim(BW1，n） 说明：计算二进制图像中对象的周长.练习31BW1 = imread(circbw。tif);BW2 = bwperim(BW1,8);imshow（BW1）figure， imshow(BW2）语法：BW2 = bwselect（BW1，c，r，n) 说明：在二进制图像中选择对象练习32BW1 = imread（text。tif);c = 16 90 144；r = 85 197 247;BW2 = bwselect(BW1,c，r,4）；imshow（BW1)figure, imshow（BW2)72 其他参数测量语法:b = std2(A)说明:计算矩阵元素的标准偏移。语法:B = mean2（A，m n）说明：计算矩阵元素的平均值附录二：MATLAB图像处理函数1applylut功能：在二进制图像中利用lookup 表进行边沿操作。语法：A = applylut（BW,lut)2bestblk功能:确定进行块操作的块大小。语法：siz = bestblk(m n，k);mb，nb = bestblk（m n,k）3blkproc功能：实现图像的显式块操作。语法：B = blkproc（A，m n，fun)；B = blkproc(A，m n,fun,P1，P2，.。)；B = blkproc（A，m n,mborder nborder，fun，。.)；B = blkproc（A，indexed，。.。)4brighten功能：增加或降低颜色映像表的亮度。语法：brighten（beta）；newmap = brighten（beta）；newmap = brighten(map，beta)；brighten（fig,beta）5bwarea功能：计算二进制图像对象的面积。语法：total = bwarea（BW）6bweuler功能：计算二进制图像的欧拉数。语法：eul = bweuler（BW，n）7bwfill功能：填充二进制图像的背景色。语法：BW2 = bwfill(BW1,c，r，n)；BW2 = bwfill（BW1，n)；BW2,idx = bwfill(.）BW2 = bwfill（x，y，BW1,xi,yi，n)；x，y，BW2,idx，xi，yi = bwfill(。.）BW2 = bwfill(BW1,holes，n）；BW2，idx = bwfill(BW1，holes,n)8bwlabel功能：标注二进制图像中已连接的部分。语法:L = bwlabel（BW,n)；L,num = bwlabel(BW，n）9bwmorph功能：提取二进制图像的轮廓.语法：BW2 = bwmorph（BW1,operation)；BW2 = bwmorph(BW1，operation,n)10bwperim功能：计算二进制图像中对象的周长.语法：BW2 = bwperim（BW1，n）11bwselect功能：在二进制图像中选择对象。语法：BW2 = bwselect（BW1,c，r，n）；BW2 = bwselect（BW1，n）BW2，idx = bwselect(.。)12cm