数字图像处理.docx
《数字图像处理.docx》由会员分享，可在线阅读，更多相关《数字图像处理.docx（24页珍藏版）》请在冰点文库上搜索。
数字图像处理
数字图像处理实验报告
学院：
专业：
班号：
学号：
学生姓名：
指导老师：
一、实验目的
锻炼编程能力，编写程序完成所给题目内容。
二、实验环境
Matlab2010b
三、实验内容、相关代码及注释、调试过程及总结
1.编程实现打开一个灰度图像，能够显示图像。
（1）解题思路：
图像读取只需要一个读取指令imread就能完成，实验中我所使用的1.jpg是彩色图像，所以需要用rgb2gray进行变换
（2）代码及注释：
clearall;
clc;
Pic=imread（'1.jpg'）;%彩色图像读取
GrayPic=rgb2gray（Pic）;%将彩色图像转换成灰度图像
subplot（1,2,1）,imshow（Pic）;%彩色图像显示
title（'彩色图像'）;
subplot（1,2,2）,imshow（GrayPic）;%灰度图像显示
title（'灰度图像'）;
（3）调试过程及总结
这仅仅是开头第一个实验，很简单很容易。
2.用灰度图形式表示彩色图象，分别显示1幅彩色图象的R，G，B分量（每个分量用8bit表示），和这幅彩色图象的H，S，I分量（每个分量也各用8bit表示）。
（1）解题思路：
首先将彩色图像用imread进行读取，再利用“x1=rgb（:
:
1）;x2=rgb（:
:
2）;x3=rgb（:
:
3）;”显示RGB分量；要显示图像的HIS分量和变换图像可以根据公式编写程序完成。
（2）代码及注释：
clearall;
clc;
rgb=imread（'1.jpg'）;
subplot（2,4,1）,imshow（rgb）;
title（'原图像'）;
%抽取图像分量
x1=rgb（:
:
1）;
x2=rgb（:
:
2）;
x3=rgb（:
:
3）;%R,G,B三个分量对某个分量或灰度图像矩阵x做傅里叶变换
subplot（2,4,2）,imshow（x1）;
title（'R分量图像'）;
subplot（2,4,3）,imshow（x2）;
title（'G分量图像'）;
subplot（2,4,4）,imshow（x3）;
title（'B分量图像'）;
%hsi=rgb2hsi（rgb）把一幅RGB图像转换为HSI图像，
%输入图像是一个彩色像素的M×N×3的数组，
%其中每一个彩色像素都在特定空间位置的彩色图像中对应红、绿、蓝三个分量。
%假如所有的RGB分量是均衡的，那么HSI转换就是未定义的。
%输入图像可能是double（取值范围是[0,1]），uint8或uint16。
%
%输出HSI图像是double，
%其中hsi（:
:
1）是色度分量，它的范围是除以2*pi后的[0,1]；
%hsi（:
:
2）是饱和度分量，范围是[0,1]；
%hsi（:
:
3）是亮度分量，范围是[0,1]。
rgb=im2double（rgb）;
r=rgb（:
:
1）;
g=rgb（:
:
2）;
b=rgb（:
:
3）;
%执行转换方程
num=0.5*（（r-g）+（r-b））;
den=sqrt（（r-g）.^2+（r-b）.*（g-b））;
theta=acos（num./（den+eps））;%防止除数为0
H=theta;
H（b>g）=2*pi-H（b>g）;
H=H/（2*pi）;
num=min（min（r,g）,b）;
den=r+g+b;
den（den==0）=eps;%防止除数为0
S=1-3.*num./den;
H（S==0）=0;
I=（r+g+b）/3;
%将3个分量联合成为一个HSI图像
hsi=cat（3,H,S,I）;
subplot（2,4,5）,imshow（hsi）;
title（'转换为HSI分量图像'）;
subplot（2,4,6）,imshow（H）;
title（'H分量图像'）;
subplot（2,4,7）,imshow（S）;
title（'S分量图像'）;
subplot（2,4,8）,imshow（I）;
title（'I分量图像'）;
（3）调试过程及总结
彩色图像的像素矩阵为一个三维矩阵，RGB分量只需按需读出各个分量就行，HSI分量的提取没有现成的函数来实现，只能按照公式进行自己别写程序来实现。
其实在本题编写的程序当中计算HIS分量部分可以编写成函数，方便以后使用。
3.编程实现图像傅立叶变换
（1）解题思路：
图像的傅里叶变换只需利用函数fft2就可完成；利用函数fftshift可以实现图像频谱的搬移，利用函数ifft2可以实现傅里叶反变换显示出变换后的图像。
（2）代码及注释：
clearall;
closeall;
clc;
%读取灰度图像
I=imread（'1.jpg'）;
I=rgb2gray（I）;
F=fft2（I）;%傅里叶变换
F1=real（F）;
Fshift=fftshift（F）;%频谱搬移
F2=real（Fshift）;
C=ifft2（Fshift）;%傅立叶反变换
subplot（2,2,1）,imshow（I）;
title（'原始图'）;%画图像原始图
subplot（2,2,2）,imshow（log（abs（F1））,[]）;
title（'直接变换频率谱'）;%画直接变换频谱图
subplot（2,2,3）,imshow（log（abs（F2））,[]）;
title（'处理后频谱'）;%平移处理后频谱图
subplot（2,2,4）,imshow（abs（C）,[]）;
title（'反变换后图像'）;%显示反傅立叶变换后的图像
（3）调试过程及总结
图像的傅里叶变换用函数fft就可以实现，利用fftshift实现频谱搬移可以更好的显示变换后的图像。
4.用直接灰度变换改变图像（求反，增强对比度，动态范围压缩，灰度切分）
（1）解题思路：
图像求反可以按照原理编写程序完成，函数imadjust可完成对比度的增强，动态范围压缩可使用im2uint8（mat2gray（log（1+double（K））））来完成，为了实现灰度切分可使用函数roicolor来完成。
（2）代码及注释：
clearall;
closeall;
clc;
%读取灰度图像
I=imread（'1.jpg'）;
K=rgb2gray（I）;
%图像求反
%x1=y1时斜率k=1的效果：
x1=256;
y1=256;
k=y1/x1;
[m,n]=size（K）;
J=double（K）;
fori=1:
m
forj=1:
n
x=J（i,j）;
y（i,j）=0;
if（x>=0）&（x<=x1）
y（i,j）=y1-k*x;
else
y（i,j）=0;
end
end
end
%对比度增强
A=imadjust（K,[01],[0.250.75]）;
%动态范围压缩
B=im2uint8（mat2gray（log（1+double（K））））;
%灰度切分
J=roicolor（K,20,200）;
subplot（2,3,1）,imshow（I）;
title（'原始图像'）;
subplot（2,3,2）,imshow（K）;
title（'灰度图像'）;
subplot（2,3,3）,imshow（mat2gray（y））;
title（'求反后图像'）;
subplot（2,3,4）,imshow（A）;
title（'对比度增强后图像'）;
subplot（2,3,5）,imshow（B）;
title（'动态压缩后图像'）;
subplot（2,3,6）,imshow（J）;
title（'灰度切分后图像'）;
（3）调试过程及总结
图像求反的时候使用了for循环，后面的对比度增强、动态压缩、灰度切分就只调用了matlab的函数就可以实现实验目的。
开始的时候还出现了问题：
经过检查发现是未将数据类型转换为double型引起的，后来得到顺利解决。
5.做出直方图，进行直方图均衡化。
（1）解题思路：
直方图反映的是图像的灰度分布，反映了灰度级和相应灰度级出现的概率之间的关系。
（2）代码及注释：
clearall;
closeall;
clc;
%一，图像的预处理，读入彩色图像将其灰度化
I=imread（'1.jpg'）;%读入JPG彩色图像文件
J=rgb2gray（I）;%灰度化后的数据存入数组
%二，绘制直方图
[m,n]=size（J）;%测量图像尺寸参数
GP=zeros（1,256）;%预创建存放灰度出现概率的向量
fork=0:
255
GP（k+1）=length（find（J==k））/（m*n）;%计算每级灰度出现的概率，将其存入GP中相应位置
end
%三，直方图均衡化
S1=zeros（1,256）;
fori=1:
256
forj=1:
i
S1（i）=GP（j）+S1（i）;%计算Sk
end
end
S2=round（（S1*256）+0.5）;%将Sk归到相近级的灰度
fori=1:
256
GPeq（i）=sum（GP（find（S2==i）））;%计算现有每个灰度级出现的概率
end
%四，图像均衡化
K=J;
fori=0:
255
K（find（J==i））=S2（i+1）;%将各个像素归一化后的灰度值赋给这个像素
end
subplot（2,2,1）,imshow（J）;
title（'灰度图像'）;
subplot（2,2,2）;
bar（0:
255,GP,'g'）%绘制直方图
title（'原图像直方图'）
xlabel（'灰度值'）
ylabel（'出现概率'）;
subplot（2,2,3）;
bar（0:
255,GPeq,'b'）%显示均衡化后的直方图
title（'均衡化后的直方图'）
xlabel（'灰度值'）
ylabel（'出现概率'）;
subplot（2,2,4）,imshow（K）;%显示均衡化后的图像
title（'均衡化后图像'）
（3）调试过程及总结
本题是在网上找到的例程加以修改实现的，后来发现matlab里面也有相应的函数histeq可以实现。
直方图的均衡化可以达到增强图像的效果。
6.用空域和频域滤波器，对图像进行平滑，并比较两种方法得到图像质量。
（1）解题思路：
图像的平滑处理用于减少噪声，但是得以付出一定细节的模糊为代价。
空域滤波可以用中值滤波和均值滤波来实现，频域常见的方法是采用巴特沃斯滤波低通滤波器来实现
（2）代码及注释：
%空域滤波
clearall;
clc;
closeall;
I=imread（'1.jpg'）;%读取图像
F=rgb2gray（I）;%转化成灰度图
J=imnoise（F,'salt&pepper',0.02）;
subplot（2,3,1）,imshow（F）;%显示灰度图
title（'灰度图像'）;
subplot（2,3,2）,imshow（J）;%显示添加椒盐噪声图
title（'添加椒盐噪声图像'）;
k1=medfilt2（J）;%3*3模板中值滤波
k2=medfilt2（J,[5,5]）;%5*5模板中值滤波
k3=medfilt2（J,[7,7]）;%7*7模板中值滤波
k4=medfilt2（J,[9,9]）;%9*9模板中值滤波
subplot（2,3,3）,imshow（k1）;
title（'3*3模板中值滤波图'）;
subplot（2,3,4）,imshow（k2）;
title（'5*5模板中值滤波图'）;
subplot（2,3,5）,imshow（k3）;
title（'7*7模板中值滤波图'）;
subplot（2,3,6）,imshow（k4）;
title（'9*9模板中值滤波图'）;
%频域滤波
clearall;
clc;
closeall;
F=imread（'1.jpg'）;%读取图像
subplot（2,2,1）,imshow（F）;
title（'原图像'）;
I=rgb2gray（F）;%转化成灰度图
subplot（2,2,2）,imshow（I）;
title（'灰度图像'）;
J=imnoise（I,'salt&pepper',0.02）;
subplot（2,2,3）,imshow（J）;
title（'加椒盐噪声后图像'）;
f=double（J）;%数据类型转换，MATLAB不支持图像的无符号整型的计算
g=fft2（f）;%傅立叶变换
g=fftshift（g）;%转换数据矩阵
[M,N]=size（g）;
nn=2;%二阶巴特沃斯（Butterworth）低通滤波器
d0=50;
m=fix（M/2）;n=fix（N/2）;
fori=1:
M
forj=1:
N
d=sqrt（（i-m）^2+（j-n）^2）;
h=1/（1+0.414*（d/d0）^（2*nn））;%计算低通滤波器传递函数
result（i,j）=h*g（i,j）;
end
end
result=ifftshift（result）;
J1=ifft2（result）;
J2=uint8（real（J1））;
subplot（2,2,4）,imshow（J2）;
title（'滤波后图像'）;%显示滤波处理后的图像
（3）调试过程及总结
调试的时候出现了错误提示
经过检查发现问题出在做实验我一直使用的是彩色图像，没有将其转为灰度图像就开始了本题的处理，函数medfilt2是不能处理彩色图像的。
7.用空域和频域滤波器，对图像进行锐化，并比较两种方法得到图像质量。
（1）解题思路：
图像锐化就是补偿图像的轮廓，增强图像的边缘及灰度跳变的部分，使图像变得清晰。
空域滤波仍然采用matlab中自带的函数=fspecial和模板Sobel实现，频域滤波采用巴特沃斯高通滤波器。
（2）代码及注释：
clearall;
clc;
closeall;
%空域滤波
F=imread（'1.jpg'）;%读取图像
I=rgb2gray（F）;%转化成灰度图
subplot（2,2,1）,imshow（I）;
title（'灰度图像'）;
H=fspecial（'Sobel'）;
H=H';%Sobel垂直模板
TH=filter2（H,I）;
subplot（2,2,2）,imshow（TH,[]）;
title（'空域水平模板'）;
H=fspecial（'Sobel'）;
H=H';%Sobel水平模板
TH=filter2（H,I）;
subplot（2,2,3）,imshow（TH,[]）;
title（'空域垂直模板'）;
%频域滤波
f=double（I）;%数据类型转换
g=fft2（f）;%傅立叶变换
g=fftshift（g）;%转换数据矩阵
[M,N]=size（g）;
nn=2;%二阶巴特沃斯（Butterworth）高通滤波器
d0=5;
m=fix（M/2）;
n=fix（N/2）;
fori=1:
M
forj=1:
N
d=sqrt（（i-m）^2+（j-n）^2）;
if（d==0）
h=0;
else
h=1/（1+0.414*（d0/d）^（2*nn））;%计算传递函数
end
result（i,j）=h*g（i,j）;
end
end
result=ifftshift（result）;
J=ifft2（result）;
J1=uint8（real（J））;
subplot（2,2,4）,imshow（J1）;%滤波后图像显示
title（'频域巴特沃斯高通滤波后结果'）;
（3）调试过程及总结
结果显示，都实现了图像的锐化处理。
本题和第6题的方式正好相反。
8选用两种压缩算法，对图像进行压缩。
并比较这两种算法。
（1）解题思路：
图像压缩有很多方式，比如：
行程长度编码、熵编码、DCT编码、小波变换编码、JEPG等。
本题采用比较简单的DCT编码和小波变换编码。
（2）代码及注释：
%DCT压缩
clearall;
clc;
closeall;
F=imread（'1.jpg'）;%读取图像
I=rgb2gray（F）;%转化成灰度图
I=im2double（I）;%数据类型转换
T=dctmtx（8）;%计算二维离散DCT矩阵
dct=@（x）T*x*T';%设置函数句柄
B=blkproc（I,[88],dct）;%图像块处理
mask=[11110000%掩膜
11100000
11000000
10000000
00000000
00000000
00000000
00000000];
B2=blkproc（B,[88],@（x）mask.*x）;%图像块处理
invdct=@（x）T'*x*T;%设置函数句柄
I2=blkproc（B2,[88],invdct）;%图像块处理
subplot（1,3,1）,imshow（F）;%显示原始图像和压缩重构图像
title（'原图像'）;
subplot（1,3,2）,imshow（I）;
title（'灰度图像'）;
subplot（1,3,3）,imshow（I2）;
title（'DTC压缩重构图像'）;
%小波压缩
clc;
clearall;
closeall;%清理工作空间
clear;
X=imread（'1.jpg'）;
X=rgb2gray（X）;
subplot（221）;imshow（X）;
title（'原始图像'）;
%对图像用小波进行层小波分解
[c,s]=wavedec2（X,2,'bior3.7'）;
%提取小波分解结构中的一层的低频系数和高频系数
cal=appcoef2（c,s,'bior3.7',1）;
ch1=detcoef2（'h',c,s,1）;%水平方向
cv1=detcoef2（'v',c,s,1）;%垂直方向
cd1=detcoef2（'d',c,s,1）;%斜线方向
%各频率成份重构
a1=wrcoef2（'a',c,s,'bior3.7',1）;
h1=wrcoef2（'h',c,s,'bior3.7',1）;
v1=wrcoef2（'v',c,s,'bior3.7',1）;
d1=wrcoef2（'d',c,s,'bior3.7',1）;
c1=[a1,h1;v1,d1];
subplot（222）,imshow（c1,[]）;
title（'分解后低频和高频信息'）;
%进行图像压缩
%保留小波分解第一层低频信息
%首先对第一层信息进行量化编码
ca1=appcoef2（c,s,'bior3.7',1）;
ca1=wcodemat（ca1,440,'mat',0）;
%改变图像高度并显示
ca1=0.5*ca1;
subplot（223）;imshow（cal,[]）;
title（'第一次压缩图像'）;
%保留小波分解第二层低频信息进行压缩
ca2=appcoef2（c,s,'bior3.7',2）;
%首先对第二层信息进行量化编码
ca2=wcodemat（ca2,440,'mat',0）;
%改变图像高度并显示
ca2=0.25*ca2;
subplot（224）;imshow（ca2,[]）;
title（'第二次压缩图像'）;
（3）调试过程及总结
实验使用matlab自带函数dctmtx和wavedec2便可完成。
在进行图像块处理的时候提示blkproc可以使用blockproc代替，于是出现了如下错误提示
查资料才发现blockproc是blkproc进一步优化而来的，按理说blockproc是可替代blkproc来使用的，可是出现了错误提示，这点是我一直没弄明白的地方。
9对图像进行边缘检测，显示边缘检测后的结果图像。
（1）解题思路：
图像的边缘是由灰度不连续性所反映的。
经典的边缘提取方法是考察图像的每个像素在某个邻域内灰度的变化，也可以是描述为亮度的变化，利用边缘邻近一阶或二阶方向导数变化规律，用简单的方法检测边缘。
（2）代码及注释：
clearall;
closeall;
clc;
%读取灰度图像
I=imread（'1.jpg'）;
I=rgb2gray（I）;
subplot（1,2,1）,imshow（I）;
title（'灰度图像'）;
[G,t]=edge（I,'log'）;
subplot（1,2,2）,imshow（G）;
title（'LOG算子检测结果'）
（3）调试过程及总结
边沿检测的算子很多：
Sobel、Rober、Prewitt、Laplacian、Canny、LOG，选用LOG算子，基本满足题目要求。
10.编写一副灰度图像的DCT变换，walsh变换，以及小波变换的结果，分别显示原始图像与变换后的图像。
（1）解题思路：
DCT变换的全称是离散余弦变换（DiscreteCosineTransform），是指将一组光强数据转换成频率数据，以便得知强度变化的情形。
Walsh函数展开式有三种：
walsh序的walsh函数、佩利序的walsh函数、hadamard序的walsh函数
小波变换是空间（时间）和频率的局部变换，能有效地从信号中提取信息。
（2）代码及注释：
%二维DCT变换
clearall;
closeall;
clc;
%读取灰度图像
I=imread（'1.jpg'）;
I=rgb2gray（I）;
subplot（2,2,1）,imshow（I）;
title（'灰度图像'）;
J=dct2（I）;
%计算二维DCT变换
subplot（2,2,2）,imshow（log（abs（J））,[]）
%图像大部分能量集中在
升级会员 