本文共 11346 字,大约阅读时间需要 37 分钟。
数字图像处理
实 验 指 导 书
东北林业大学信息与计算机工程学院
计算机科学与技术专业
目 录
1 实验目的与要求………………………………………………………………… 2
2 实验环境………………………………………………………………………… 2
3实验一般步骤……………………………………………………………………2
4 实验时数………………………………………………………………………………3
5 实验内容和具体要求…………………………………………………………………3
实验一 图像打开、保存与显示…………………………………………………………4
实验二 图像灰度直方图统计…………………………………………………………6
实验三 图像轮廓提取与边缘检测……………………………………………………8
实验四 二维图像几何变换……………………………………………………………11
补充实验 均值滤波与中值滤波………………………………………………………14
1、实验目的与要求
本课程是信息管理与信息系统专业的选修课程。通过教学使学生了解数字图像处理技术的发展概况,使学生掌握图像的数字化、图像的变换、图像的压缩和编码、图像的增强、图像恢复和重建等基本概念和相关方法,为在这一领域中进行深入学习和研究打下必要的基础。上机操作是本课程必不可少的实践环节,主要目的是锻炼和培养学生实际操作技能和综合利用所学过的知识解决实际问题的能力。
(1)加深对上课内容的理解。因为缺少实际操作会使同学难以记住,通过多次上机,就能自然地、熟练地巩固数字图像处理的相关技术和方法。
(2)熟悉所用计算机系统和软件的操作方法,也就是了解和熟悉Matlab环境。掌握图像处理工具箱的基本函数和命令。
(3)学会上机调试程序,善于发现错误,并且能很快地排除错误,使程序能正确运行。计算机技术是实践性很强的技术,要求从事这一领域的人不仅能了解和熟悉有关理论和方法,还要求自己动手实现。因此调试程序本身是程序设计课程的一个重要的内容和基本要求,应给予充分的重视。
2、实验环境
计算机、Windows XP操作系统、MATLAB7.0软件
3、实验一般步骤
(1)预先准备好上机所需的程序。手编程序应书写整齐,并经人工检查无误后才能上机,以提高上机效率。对程序中自己有疑问的地方,应做出记号,以便在上机时给予注意。
(2)上机输入和调试程序。应该一人一组,独立上机。上机过程中出现的问题,除了是系统的问题以外,一般应自己独立处理,不要轻易举手问教师。尤其对“出错信息”,应善于自己分析判断。这是学习调试程序的良好机会。在程序调试通过后,打印输出程序清单,在运行时要注意在输入不同图像文件时所得到的不同结果。此时应运行几次,分别检查在不同情况下程序是否正确。
(3)上机结束后,应及时整理出实验报告,实验报告包括以下内容:
① 实验目的
② 实验环境
③ 实验内容
④ 数据处理与讨论分析
4、实验时数
总实验时数8学时。
5、实验内容和具体要求
实验一 图像打开、保存与显示
一、实验目的:
掌握数字图像的基本类型及其表示。熟悉Matlab软件环境,了解Matlab中对图像数据的读入、显示和输出等操作,实现图像文件的打开、显示与保存功能。
二、实验环境:
计算机、Windows XP操作系统,Matlab7.0
三、实验内容:
1、运用Matlab图像处理工具箱中的imread函数分别读入灰度图像、二值图像、RGB图像和索引图像,并观察相应的图像矩阵,并运用imshow函数显示相应图像。
%读入灰度图像pout,并显示
I_huidu=imread('pout.tif');
figure(1), imshow(I_huidu) , title('灰度图像');
%读入二值图像circles,并显示,注意查看图像矩阵(uint8,0表示黑,255表示白)
I_erzhi=imread('circles.png');
figure(2), imshow(I_erzhi) , title('二值图像');
%读入二值图像blobs,并显示,注意查看图像矩阵(logical,0表示黑,1表示白)
I_erzhi=imread('blobs.png');
figure(2), imshow(I_erzhi) , title('二值图像');
%读入RGB图像football,并显示
I_RGB=imread('football.jpg');
figure(3), imshow(I_RGB) , title('RGB图像');
%读入索引图像trees,并显示
[I_suoyin, colormap]=imread('trees.tif');
figure(4), imshow(I_suoyin, colormap) , title('索引图像');
2、对一个RGB彩色图像,分别抽取其R、G、B三个分量层,并显示各层图像。
A=imread('football.jpg'); %也可读入RGB图像peppers.png
R=A(:, :, 1);
G=A(:, :, 2);
B=A(:, :, 3);
subplot(2, 2, 1), imshow(A), title('原始图像');
subplot(2, 2, 2), imshow(R), title('R层图像');
subplot(2, 2, 3), imshow(G), title('G层图像');
subplot(2, 2, 4), imshow(B), title('B层图像');
3、向灰度图像pout中分别加入高斯噪声和椒盐噪声,显示并保存带有噪声的图像。
I=imread('pout.tif ');
subplot(1, 3, 1), imshow(I), title('原始图像');
G=imnoise(I, 'gaussian');
subplot(1, 3, 2), imshow(G), title('高斯噪声图像');
imwrite(G, 'Gpout.tif');
J=imnoise(I, 'salt & pepper', 0.05);
subplot(1, 3, 3), imshow(J), title('椒盐噪声图像');
imwrite(J, 'Jpout.tif');
注意:
(1) 为方便用户学习和使用Matlab图像工具箱,Matlab为用户提供了丰富的图像资源,所使用的图像文件都保存在Matlab安装目录下的\toolbox\images\imdemos子目录下。
(2) 将每次实验中需要使用的图像文件,最好在实验开始前事先拷贝到Matlab安装目录下的\work子目录下。(因为Matlab默认处理当前工作目录下的图像文件)
(3) 实验过程中,学会使用help命令。
实验二 图像灰度直方图统计
一、实验目的:
理解并掌握灰度直方图的概念、计算灰度直方图的方法以及如何应用直方图均衡化来增强图像对比度。
二、实验环境:
计算机、Windows XP操作系统,Matlab7.0
三、实验内容:
以灰度图像pout为例,运用Matlab编程实现灰度直方图的统计以及直方图均衡化处理过程:
(1)计算并绘制原始pout图像的灰度直方图;
(2)根据离散累计分布函数,对原始灰度直方图进行均衡化处理,绘制均衡化后的灰度直方图;
(3)生成均衡化处理后的新图像,显示并保存。
(4)比较原始pout图像和新图像的对比度。
clc;
clear;
grayimage=imread('pout.tif'); %读取一幅256级灰度图像,作为原始图像
subplot(2,2,1); imshow(grayimage); title('原始图像'); %显示原始图像
[m,n]=size(grayimage); %获取图像大小
%以下计算原始图像中各灰度级出现的概率
gp=zeros(1,256); %初始化一维行向量gp
for i=1:256
gp(i)=length(find(grayimage==(i-1)))/(m*n); %计算每个灰度级出现的概率
end
subplot(2,2,2); bar(0:255, gp); %显示原始图像的灰度直方图
title('原始图像的灰度直方图');
xlabel('灰度值'); ylabel('出现概率');
%以下进行直方图均衡化处理,根据离散累计分布函数,求出所有新的灰度级
S1=zeros(1,256);
S2=zeros(1,256);
temp=0;
for i=1:256
temp=temp+gp(i);
S1(i)=temp;
end
S2=round(S1*255); %将Sk归到相近的灰度级上
%以下计算各个新的灰度级出现的概率
newgp=zeros(1,256);
for i=1:256
newgp(i)=sum(gp(find(S2==(i-1))));
end
subplot(2,2,3); bar(0:255, newgp); %显示均衡化后的灰度直方图
title('均衡化后的灰度直方图');
xlabel('灰度值'); ylabel('出现概率');
%用新的灰度值填充各像素点,从而获得均衡化处理之后的新图像
newgrayimage=grayimage;
for i=1:256
newgrayimage(find(grayimage==(i-1)))=S2(i);
end
subplot(2,2,4); imshow(newgrayimage); %显示均衡化处理后的新图像
title('均衡化后的图像');
imwrite(newgrayimage,'newpout.tif'); %保存均衡化处理后的新图像
实验三 图像轮廓提取与边缘检测
一、实验目的:
理解并掌握对二值图像进行轮廓提取的相关算法(比如,掏空内部点法),以及用于图像边缘检测和提取的典型微分算子(梯度算子和拉普拉斯算子)。
二、实验环境:
计算机、Windows XP操作系统,Matlab7.0
二、实验内容:
1、根据掏空内部点算法,运用Matlab编程实现二值图像的轮廓提取。
%以下适用于黑色背景白色前景的二值图像轮廓提取(以二值图像circles为例)
BW=imread('circles.png'); %二值图像circles是uint8,0黑,255白
subplot(1,2,1); imshow(BW); title('二值图像');
[M, N]=size(BW); %M行,N列
BW_Buffer=BW;
for i=2: M-1
for j=2: N-1
if (BW(i, j)==255 & BW(i-1, j)==255 & BW(i+1, j)==255 & BW(i, j-1)==255 & BW(i, j+1)==255 & BW(i-1, j-1)==255 & BW(i-1, j+1)==255 & BW(i+1, j-1)==255 & BW(i+1, j+1)==255) %说明BW(i, j)是前景中的一个内部点
BW_Buffer(i, j)=0; %掏空该内部点,即将该内部点置成与背景相同灰度
end
end
end
subplot(1,2,2); imshow(BW_Buffer); title('提取轮廓');
%以下适用于白色背景黑色前景的二值图像轮廓提取(以二值图像source为例)
BW=imread('source.bmp'); %二值图像source是uint8,0黑,255白
subplot(1,2,1); imshow(BW); title('二值图像');
[M, N]=size(BW); %M行,N列
BW_Buffer=BW;
for i=2: M-1
for j=2: N-1
if (BW(i, j)==0 & BW(i-1, j)==0 & BW(i+1, j)==0 & BW(i, j-1)==0 & BW(i, j+1)==0 & BW(i-1, j-1)==0 & BW(i-1, j+1)==0 & BW(i+1, j-1)==0 & BW(i+1, j+1)==0) %说明BW(i, j)是前景中的一个内部点
BW_Buffer(i, j)=255; %掏空该内部点,即将该内部点置成与背景相同灰度
end
end
end
subplot(1,2,2); imshow(BW_Buffer); title('提取轮廓');
注意:使用掏空内部点的方法来提取二值图像的轮廓时,不能直接在原始二值图像矩阵上判断一个点掏空一个点,否则对前面像素的掏空操作会影响到对后面像素的判断结果。
解决方法:创建原始二值图像矩阵的副本(即图像矩阵BW_Buffer),在原始二值图像矩阵上执行判断操作,即依次判断每个像素点是否为前景中的内部点,如果是,则在图像矩阵BW_Buffer上执行掏空内部点的操作。
2、以灰度图像rice和cameraman为例,利用Matlab图像处理工具箱中的edge函数,分别使用Roberts 算子、Sobel算子、Prewitt 算子对其进行边缘检测。
(1)函数格式: BW = edge(I, 'method', thresh)
(2)格式说明:edge函数输入灰度图像矩阵I,输出二值图像矩阵BW;参数'method'用于指定所使用的边缘检测算子,可以是'roberts'、'sobel'、'prewitt'、'log'、'canny';参数thresh用于指定梯度门限值(也称梯度阈值),图像中梯度值大于等于门限值thresh的像素用白色(1)表示,说明这些地方对应边缘,梯度值小于门限值thresh的像素用黑色(0)表示(edge function will ignore all edges that are not stronger than thresh)。若不指定参数thresh,则edge函数会自动选择阈值。所以edge函数最终将原始灰度图像中的边缘和背景用二值图像的形式展现出来,以突出边缘的位置,达到边缘检测的目的。
(3)程序如下:
I=imread('rice.png');
subplot(2,2,1); imshow(I); title('原始图像');
[BW1,thresh1]=edge(I,'roberts'); %进行Roberts算子边缘检测并返回门限值
[BW2,thresh2]=edge(I,'sobel'); %进行Sobel算子边缘检测并返回门限值
[BW3,thresh3]=edge(I,'prewitt'); %进行Prewitt算子边缘检测并返回门限值
subplot(2,2,2); imshow(BW1); title('Roberts算子边缘检测结果');
subplot(2,2,3); imshow(BW2); title('Sobel算子边缘检测结果');
subplot(2,2,4); imshow(BW3); title('Prewitt算子边缘检测结果');
若向原始图像中加入随机噪声(比如高斯噪声),之后再对噪声图像分别运用Roberts 算子、Sobel算子、Prewitt 算子、Log算子(高斯-拉普拉斯算子)进行边缘检测,观察检测结果,试比较4种边缘检测算子的抗噪声干扰能力。
I=imread('rice.png');
subplot(2,3,1); imshow(I); title('原始图像');
G=imnoise(I, 'gaussian'); %向原始图像中加入高斯噪声
subplot(2,3,2); imshow(G); title('噪声图像');
BW1=edge(G, 'roberts'); %进行Roberts算子边缘检测
BW2=edge(G, 'sobel'); %进行Sobel算子边缘检测
BW3=edge(G, 'prewitt'); %进行Prewitt算子边缘检测
BW4=edge(G, 'log'); %进行Log算子边缘检测
subplot(2,3,3); imshow(BW1); title('Roberts算子边缘检测结果');
subplot(2,3,4); imshow(BW2); title('Sobel算子边缘检测结果');
subplot(2,3,5); imshow(BW3); title('Prewitt算子边缘检测结果');
subplot(2,3,6); imshow(BW4); title('Log算子边缘检测结果');
实验四 二维图像几何变换
一、实验目的:
理解并掌握平移、镜像、旋转和比例缩放等二维图像基本几何变换的变换原理以及灰度插值处理过程。
二、实验环境:
计算机、Windows XP操作系统,Matlab7.0
三、实验内容:
1、以灰度图像cameraman为例,通过Matlab编程实现水平镜像和垂直镜像变换。
I=imread('cameraman.tif');
subplot(2,2,1); imshow(I); title('原始图像');
[M, N]=size(I); %求出原始图像的高度(行数M)和宽度(列数N)
%以下实现水平镜像
I_shuiping=I;
for i=1: M
for j=1: N
I_shuiping(i, j) = I(i, N+1-j);
end
end
subplot(2,2,2); imshow(I_shuiping); title('水平镜像图像');
%以下实现垂直镜像
I_chuizhi=I;
for j=1: N
for i=1: M
I_chuizhi(i, j) = I(M+1-i, j);
end
end
subplot(2,2,3); imshow(I_chuizhi); title('垂直镜像图像');
思考:如果将一幅RGB图像分别进行水平镜像和垂直镜像,应该如何处理?
I=imread('peppers.png'); %读入一幅RGB彩色图像,以peppers.png为例
subplot(2,2,1); imshow(I); title('原始图像');
[M, N, L]=size(I); %求出原始图像的大小(M行,N列,L=3层)
%以下实现水平镜像
I_shuiping=I;
for k=1:L
for i=1:M
for j=1:N
I_shuiping(i, j, k) = I(i, N+1-j, k);
end
end
end
subplot(2,2,2); imshow(I_shuiping); title('水平镜像图像');
%以下实现垂直镜像
I_chuizhi=I;
for k=1:L
for j=1:N
for i=1:M
I_chuizhi(i, j, k) = I(M+1-i, j, k);
end
end
end
subplot(2,2,3); imshow(I_chuizhi); title('垂直镜像图像');
2、以灰度图像circuit为例,利用Matlab图像处理工具箱中的imresize函数对图像进行比例缩放变换。
(1)函数格式:
B=imresize(A, m, method)
B=imresize(A, [mrows ncols], method)
(2)格式说明:A为原始图像矩阵,B为比例缩放变换后的图像矩阵;method指定缩放过程中使用的插值方法,“nearest”为最邻近插值法(默认值),“bilinear”为双线性插值法;m为缩放倍数,m大于1,图像被放大,小于1,图像被缩小;mrows和ncols用于分别指定变换后图像的高度和宽度,这能克服缩放倍数m只能对高度和宽度等比例缩放的缺陷。
(3)程序如下:
I=imread('circuit.tif');
[M, N]=size(I); %求出原始图像的高度(行数M)和宽度(列数N)
F=imresize(I, 1.5); %将原始图像等比例放大为原来的1.5倍
imwrite(F, 'D:\circuit1_5.tif'); %保存放大1.5倍后的图像
[M1_5, N1_5]=size(F); %求出放大1.5倍后图像的高度和宽度
S=imresize(I, 0.5); %将原始图像等比例缩小为原来的0.5倍
imwrite(S, 'D:\circuit0_5.tif'); %保存缩小0.5倍后的图像
[M0_5, N0_5]=size(S); %求出缩小0.5倍后图像的高度和宽度
J=imresize(I, [190,400]); %将原始图像缩放为高190宽400的图像,实现不等比例缩放变换
imwrite(J, 'D:\circuit190_400.tif'); %保存不等比例缩放变换后的图像
注意:可以创建4个figure窗口,分别用于显示原始图像、等比例放大1.5倍后的图像、等比例缩小0.5倍后的图像、缩放为高190宽400的图像。但一定要使各自的figure窗口保持初始大小,不要最大化figure窗口,因为figure窗口的初始大小才表明了缩放前后图像大小的变化。这里也可以使用其它图像处理工具软件(比如ACDSee),分别打开上述程序中的原始图像以及新保存的3幅比例缩放变换后的图像,以观察它们彼此之间的变化。
3、以灰度图像cameraman为例,利用Matlab图像处理工具箱中的imrotate函数对图像进行旋转变换。
(1)函数格式:
B=imrotate(A, angle)
B=imrotate(A, angle, method)
(2)格式说明:A为原始图像矩阵,B为旋转变换后的图像矩阵;angle指定逆时针旋转角度,当angle为负值时,表示顺时针旋转;method指定旋转过程中使用的插值方法,“nearest”为最邻近插值法(默认值),“bilinear”为双线性插值法。
(3)程序如下:
I=imread('cameraman.tif');
I_45=imrotate(I, 45); %将原始图像逆时针旋转45度
imwrite(I_45, 'D:\cameraman45.tif'); %保存旋转后的图像
I_45N=imrotate(I, -45); %将原始图像顺时针旋转45度
imwrite(I_45N, 'D:\cameraman45N.tif'); %保存旋转后的图像
I_180=imrotate(I, 180); %将原始图像逆时针旋转180度
imwrite(I_180, 'D:\cameraman180.tif'); %保存旋转后的图像
注意:可以创建4个figure窗口,分别用于显示原始图像、逆时针旋转45度后的图像、顺时针旋转45度后的图像、逆时针旋转180度后的图像。但一定要使各自的figure窗口保持初始大小,不要最大化figure窗口,因为figure窗口的初始大小才表明了旋转前后图像的变化。这里也可以使用其它图像处理工具软件(比如ACDSee),分别打开上述程序中的原始图像以及新保存的3幅旋转变换后的图像,以观察它们彼此之间的变化。
补充实验 均值滤波与中值滤波
一、均值滤波
(1)avefilt.m
% 自编的均值滤波函数
% x是要进行均值滤波处理的图像矩阵
% n是均值滤波中使用的模板大小(即n×n)
% d是均值滤波处理后生成的图像矩阵
function d=avefilt(x, n)
a=ones(n); % a为大小为n×n的模板,元素全是1,源码a(1:n, 1:n) =1;
[M, N]=size(x); %图像矩阵为M行×N列
x1=double(x);
x2=x1;
% A(a:b, c:d)表示取矩阵A的第a行到第b行第c列到第d列的所有元素,构成一个新矩阵
for i=1:M-n+1
for j=1:N-n+1
c=x1(i:i+(n-1), j:j+(n-1)).*a; %取出x1中从(i, j)开始的n行n列元素与模板相乘
s=sum(sum(c)); %求c矩阵中所有元素之和
x2(i+(n-1)/2, j+(n-1)/2)=s/(n*n); %将均值赋给中心元素,未被赋值的元素取原值,即忽略边界像素点的处理
end
end
d=uint8(x2);
--------------------------------------------------------------------------------------------------------------
(2)junzhilvbo.m
clc;
clear;
original=imread('woman.bmp');
G_original=imnoise(original, 'gaussian');
subplot(1,3,1), imshow(original), title('原始图像');
subplot(1,3,2), imshow(G_original), title('加入高斯噪声后的图像');
after=avefilt(G_original, 3); %调用avefilt函数
subplot(1,3,3), imshow(after), title('均值滤波处理后的图像');
二、中值滤波
(1)midfilt.m
% 自编的中值滤波函数
% x是需要进行中值滤波处理的图像矩阵
% n是中值滤波中使用的窗口大小(即n×n)
% d是中值滤波处理后生成的图像矩阵
function d=midfilt(x, n)
[M, N]=size(x); %图像矩阵为M行×N列
x1=double(x);
x2=x1;
for i=1:M-n+1
for j=1:N-n+1
c=x1(i:i+(n-1), j:j+(n-1)); %取出x1中从(i, j)开始的n行n列元素
e=c(1,:);
for k=2:n
e=[e, c(k,:)]; %将矩阵c中的所有元素形成一个行向量
end
x2(i+(n-1)/2, j+(n-1)/2)=median(e); %将中值赋给中心元素,未被赋值的元素取原值,即忽略边界像素点的处理
end
end
d=uint8(x2);
--------------------------------------------------------------------------------------------------------------
(2)zhongzhilvbo.m
clc;
clear;
original=imread('woman.bmp');
J_original=imnoise(original, 'salt & pepper');
subplot(1,3,1), imshow(original), title('原始图像');
subplot(1,3,2), imshow(J_original), title('加入椒盐噪声后的图像');
after=midfilt(J_original, 3); %调用midfilt函数
subplot(1,3,3), imshow(after), title('中值滤波处理后的图像');
转载地址:http://fvali.baihongyu.com/