中国地质大学(武汉)
数字图像处理实验报告
姓 名: 张彪 _
班 号: 075112 _
学 号:20111002253
院 系:_机电学院
专 业:_通信工程
指导教师:_李杏梅老师
20## 年4月
第一次实验
一、
1.实验内容:
根据灰度图象得到一副彩色图像(变换函数自定),分别显示1幅彩色图象的R,G,B分量(每个分量用8 bit表示),和这幅彩色图象的H,S,I分量(每个分量也各用8 bit表示)。
2.实验目的:
进一步掌握彩色图像处理知识,重点是掌握RGB彩色模型到HIS彩色模型的转换。
3.实验原理:
HIS [Hue-Saturation-Intensity(Lightness),HSI或HSL] 颜色模型用H、S、I三参数描述颜色特性,其中H定义颜色的波长,称为色调;S表示颜色的深浅程度,称为饱和度;I表示强度或亮度
根据三基色原理,用基色光单位来表示光的量,则在RGB颜色空间,任意色光F都可以用R、G、B三色不同分量的相加混合而成:F=r[R]+g[G]+b[B]
RGB向HSI模型的转换是由一个基于笛卡尔直角坐标系的单位立方体向基于圆柱极坐标的体的转换。基本要求是将RGB中的亮度因素分离将色度分解为色调和饱和度,并用角向量表示色调。具体转换公式教科书上皆有参考。
4.实验思路:
此实验的实验思路较为简单清晰,即先将一副彩色图分离出R、G、B分量,然后利用公式转化为对应的H、S、I分量即可,中间遇到问题较少。
5.实验代码、注释:
图像的读取:
clear all;
clc;
rgb = imread('1.jpg');
subplot(2,4,1),imshow(rgb);
title('原图像');
RGB分量的提取:
% 抽取图像分量
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;
转换为HIS分量并输出图像:
% 将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分量图像');
6、实验运行结果:
小结:
此题只要清楚实验原理。完成过程保持细心,很少出现大的问题。
二、
1.实验内容:
编写一副灰度图像的DCT变换,walsh变换,以及小波变换的结果,分别显示原始图像与变换后的图像。
2.实验目的:
掌握图像处理中图像压缩编码中的相关知识,学会利用DCT、walsh及小波变换等方式对图像进行相应的变换,进一步了解想换知识与应用。
3.实验原理:
DCT变换 DCT变换的全称是离散余弦变换(Discrete Cosine Transform),是指将一组光强数据转换成频率数据,以便得知强度变化的情形。若对高频的数据做些修饰,再转回原来形式的数据时,显然与原始数据有些差异,但是人类的眼睛却是不容易辨认出来。
沃尔什变换主要用于图像变换,属于正交变换。这种变换压缩效率低,所以实际使用并不多。但它快速,因为计算只需加减和偶尔的右移操作。沃尔什变换的定义如下:给定一个NXN像素块Pxy(N必须是2的幂),二维WHT定义为:
小波变换是时间(空间)频率的局部化分析,它通过伸缩平移运算对信号(函数)逐步进行多尺度细化,最终达到高频处时间细分,低频处频率细分,能自动适应时频信号分析的要求,从而可聚焦到信号的任意细节,解决了Fourier变换的困难问题。
4. 实验思路:
DCT、小波变换在了解了原理之后,即可按步骤实现变换,过程很简单。而walsh变换最初思路是先读取图像产生一个矩阵然后作沃尔什-哈达玛变换,最后得到变换。后来自己编写的程序出现问题,然后和同学一起讨论用for循环实现变换,过程更清晰简单,给老师验收时变采用的此思路下的结果。
5.实验代码、注释:
DCT变换:
%二维DCT变换
clear all;
close all;
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)),[])
%图像大部分能量集中在上左角处
%colormap(jet(64)), colorbar
title('DCT变换后图像');
J(abs(J) < 10) = 0; %上述DCT 变换结果中绝对值小于10 的系数舍弃
%colormap(jet(64)), colorbar
%用idct2重构图像
K = idct2(J)/255;
subplot(2,2,3),imshow(K);
title('重构图像');
Walsh变换:
clear all;
close all;
a=imread('1.jpg');
A=rgb2gray(a);
subplot(121);
imshow(A);
title('原图')
A=double(A);
I=hadamard(256);
[m,n]=size(I);
count=0;
e(1:m)=0
y=zeros(m,n)
for a=1:m
for b=1:n
if(b==n)
break;
else
if(abs( I(a,b)-I(a,b+1))==2)
count=count+1
end
end
end
e(a)=count;
y(count+1,:)=I(a,:);
count=0;
end
l=y*A*y/(256*256);
subplot(122);
imshow(l);title('walsh变换');
小波变换
% 小波变换
clear all;
close all;
clc;
%读取灰度图像
I=imread('1.jpg');
I = rgb2gray(I);
sI=size(I);
[cA1,cH1,cV1,cD1] = dwt2(I,'db2');
X = idwt2(cA1,cH1,cV1,cD1,'db4',sI);
figure()
subplot(2,2,1),imshow(cA1);
title('近似细节系数');
subplot(2,2,2),imshow(cH1);
title('水平细节系数');
subplot(2,2,3),imshow(cV1);
title('垂直细节系数');
subplot(2,2,4),imshow(cD1);
title('对角细节系数');
X=uint8(X);
figure()
imshow(X);
title('重构后图像')
6、实验运行结果:
(1) DCT变换运行结果:
(2)walsh变换结果:
(3)小波变换结果:
7、问题小结:
在实验过程中,也出现了一些问题,但是最后经查阅资料,同学讨论,大都得以解决,在验收过程中,老师问了许多问题,有一个当时自己并不是很清楚,就是DCT变换中的关于系数的舍弃,课下的时候自己查阅资料和老师讨论了一下,认识更清晰了,当然老师验收时说原图和重构图不一样,我还以为老师看错了,其实确实是不一样,因为重构是在舍弃部分系数情况下的结果,所以结果不同,但一般差别较为微小,老师如此细心认真,一丝不苟,严谨的态度,着实让学生受益匪浅。关于其中处理的应用,和老师讨论了一些,也发现自己许多地方需要改进,向老师多多学习,感谢老师的悉心指导。
三、
1、实验内容:
编程实现图像傅立叶变换,并且显示图像的频谱图,编程实现傅里叶变换的,分离性,平移,旋转,卷积性质。
2、实验目的:
进一步了解傅利叶变换基本知识;
学会用matlab实现傅里叶基本性质变换。
3、实验代码及注释:
clear all;
close all;
clear;clc;
I=imread('1.jpg');
I=double(I);a=I(:,:,1);
%旋转性质
aa=imrotate(a,45,'bilinear','loose');
b=zeros(313);
b(47:266,47:266)=a;
afft=fftshift(abs(fft2(aa)));
bfft=fftshift(abs(fft2(b)));
figure(2)
subplot(222)
imshow(aa,[]);title('原图')
subplot(221)
imshow(b,[]);title('旋转45度')
subplot(224)
imshow(log(afft),[]);title('原频域')
subplot(223)
imshow(log(bfft),[]);title('旋转45度频域')
%离散性
ax=fft(a);
ay=fft(ax');
ay=ay';
ay=fftshift(abs(ay));
af=fftshift(abs(fft2(a)));
figure(3)
subplot(121)
imshow(log(af),[]);title('原频谱');
subplot(122)
imshow(log(ay),[]);title('离散性频谱');
%平移性
k1=fft2(a);
for m=1:220
for n=1:220
a(m,n)=a(m,n)*(-1)^(m+n);
end
end
k2=fft2(a);
figure(4)
subplot(121)
imshow(log(abs(k1)),[]);title('原图傅立叶变换');
subplot(122)
imshow(log(abs(k2)),[]);title('平移性')
%卷积性质
c=a(1:128,1:128);
d=a(21:148,21:148);
m=conv2(c,d);
c(383,383)=0;
d(383,383)=0;
m1=ifft2((fft2(c).*fft2(d)));
m1=real(m1);
m1=m1(1:255,1:255);
figure(1)
subplot(121)
imshow(log(m),[]); title('时域卷积')
subplot(122)
imshow(log(m1),[]); title('频域相乘后反变换')
4、实现结果:
5、实验小结:
这个实验在运行过程中卷积这一块出现了错误,后来发现是因为连续做了两次变换,所以最后反变换之后两个图像并不同,后来稍作修改便得出正确图像了。和老师也交流了许多,确实收获了很多。感谢老师的悉心指导。
第二次实验
四、
1、实验内容:
用直接灰度变换改变图像(求反,增强对比度,动态范围压缩(对数变换),灰度切分),显示一副灰度图像的8个位面图。
2、实验目的:
进一步了解图像增强的应用及其分类的基本知识;
学会用matlab处理图像,实现图像增强效果。
3.实验原理
(1) 求反对图象求反是将原图灰度值翻转。
(2) 增强对比度:实际是增强原图各部分之间的反差(灰度差别)。
(3) 动态范围压缩:该方法的目标与增强对比度相反。灰度取值小的区间,灰度范围增加,取值大的区间,灰度范围减小。
(4) 灰度切分:它的目的和增强对比度相似,要将灰度值范围变得比较突出。关心的范围取较高的值。
(5) 位面图:每个像素的灰度等级是256,这可用8位来表示,假设图像是由8个1位面组成,范围从0到位面7.位面0包含图像像素最低位,位面7包含像素最高位。
4、验代码及注释:
(1)、(2)、(3)、(4)代码:
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);
for i =1:m
for j=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,[0 1],[0.25 0.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('灰度切分后图像');
(5)位图实现代码:
clear all
a=imread('1.jpg');
subplot(3,3,1);
imshow(a);
title('原图');
a=rgb2gray(a);
a=double(a);
[M,N] = size(a);
g=a;
for i=1:M
for j=1:N
if bitand(g(i,j),bitshift(1,7))>0
g(i,j)=255;
else g(i,j)=0;
end
end
end
subplot(3,3,2);
imshow(g);title('第7位图);
h=a;
for i=1:M
for j=1:N
if bitand(h(i,j),bitshift(1,6))>0
h(i,j)=255;
else h(i,j)=0;
end
end
end
subplot(3,3,3);imshow(h);title('第6位图);
k=a;
for i=1:M
for j=1:N
if bitand(k(i,j),bitshift(1,5))>0
k(i,j)=255;
else k(i,j)=0;
end
end
end
subplot(3,3,4);imshow(k);title('第5位图);
z=a;
for i=1:M
for j=1:N
if bitand(z(i,j),bitshift(1,4))>0
z(i,j)=255;
else z(i,j)=0;
end
end
end
subplot(3,3,5);imshow(z);title('第4位图);
c=a;
for i=1:M
for j=1:N
if bitand(c(i,j),bitshift(1,3))>0
c(i,j)=255;
else c(i,j)=0;
end
end
end
subplot(3,3,6);imshow(c);title('第3 位图);
v=a;
for i=1:M
for j=1:N
if bitand(v(i,j),bitshift(1,2))>0
v(i,j)=255;
else v(i,j)=0;
end
end
end
subplot(3,3,7);imshow(v);title('第2位图);
n=a;
for i=1:M
for j=1:N
if bitand(n(i,j),bitshift(1,1))>0
n(i,j)=255;
else n(i,j)=0;
end
end
end
subplot(3,3,8);imshow(n);title('第1位图);
s=a;
for i=1:M
for j=1:N
if bitand(s(i,j),bitshift(1,0))>0
s(i,j)=255;
else s(i,j)=0;
end
end
end
subplot(3,3,9);
imshow(s);title('第0图');
4、实验运行结果:
五、
1、实验内容:
编程实现均值,中值,最大值滤波,编程实现各个高通滤波的各种算子的边缘检测。
2、实验目的:
进一步了解图像滤波的应用及其算子边缘检测的基本知识;
学会用matlab处理图像,实现滤波以及算子边缘检测等功能。
3.实验原理
均值滤波的主要步骤为:
(1) 将模板在途中漫游,并将模板中心与途中某个象素位置重合;
(2) 将模板上系数与模板下对应象素相乘;
(3) 将所有乘积相加;
(4) 将和(模板的输出响应)赋给途中对应模板中心位置的象素。
中值滤波的主要步骤为:
(1)将模板在途中漫游,并将模板中心与途中某个象素位置重合;
(2)读取模板下各对应象素的灰度值;
(3)将这些灰度值从小到大排成1列;
(4)找出这些值里排在中间的1个;
(5)将这个中间值赋给对应模板中心位置的象素。
最大值滤波同理
4、验代码及注释:
(1)均值实现运行程序:
clear all;
%中值均值滤波实现程序
h=imread('tu.jpg'); %读入彩色图片
c=rgb2gray(h); %把彩色图片转化成灰度图片
figure;
imshow(c);
title('原始图象'); %显示原始图象
g=imnoise(c,'gaussian',0.1,0.002); %加入高斯噪声
figure;
imshow(g);
title('加入高斯噪声之后的图象');
Y2=avefilt(g,3); %调用自定义进行均值滤波,n为模板大小
figure;
imshow(Y2);
title('均值滤波之后的结果');
Y4=midfilt(g,3); %调用自定义函数进行中值滤波
figure;
imshow(Y4);
title('中值滤波之后的结果');
(2)均值自定义函数:
function d=avefilt(x,n)
%均值自定义函数
a(1:n,1:n)=1; %a即n×n模板,元素全是1
p=size(x); %输入图像是p×q的,且p>n,q>n
x1=double(x);
x2=x1;
%A(a:b,c:d)表示A矩阵的第a到b行,第c到d列的所有元素
for i=1:p(1)-n+1
for j=1:p(2)-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);
(3)中值自定义函数:
function d=midfilt(x,n)
%中值自定义函数
p=size(x); %输入图像是p×q的,且p>n,q>n
x1=double(x);
x2=x1;
for i=1:p(1)-n+1
for j=1:p(2)-n+1
c=x1(i:i+(n-1),j:j+(n-1)); %取出x1中从(i,j)开始的n行n列元素,即模板(n×n的)
e=c(1,:); %是c矩阵的第一行
for u=2:n
e=[e,c(u,:)]; %将c矩阵变为一个行矩阵
end
mm=median(e); %利用函数median取中值
x2(i+(n-1)/2,j+(n-1)/2)=mm; %将模板各元素的中值赋给模板中心位置的元素
end
end
%未被赋值的元素取原值
d=uint8(x2);
(4)最大值滤波程序:
clear all;
%最大值滤波程序
c=imread('1.jpg');
I=rgb2gray(c); %把彩色图片转化成灰度图片
J=imnoise(I,'salt & pepper',0.02);
%添加椒盐噪声
K=ordfilt2(J,9,ones(3,3));
%3x3邻域窗最大值滤波
subplot(121),imshow(J);
title('添加椒盐噪声图像')
subplot(122),imshow(K);
title('3x3邻域窗的最大值滤波图像')
(5)算子边缘检测:
clear all
Q=imread('1.jpg');% 提取图像
I=rgb2gray(Q);
BW1=edge(I,'sobel'); %用SOBEL算子进行边缘检测
BW2=edge(I,'roberts');%用Roberts算子进行边缘检测
BW3=edge(I,'prewitt'); %用prewitt算子进行边缘检测
BW4=edge(I,'log'); %用log算子进行边缘检测
BW5=edge(I,'canny'); %用canny算子进行边缘检测
subplot(2,3,1), imshow(BW1);
title('sobel 算子检测');
subplot(2,3,2), imshow(BW2);
title('Roberts 算子检测');
subplot(2,3,3), imshow(BW3);
title('prewitt 算子检测');
subplot(2,3,4), imshow(BW4);
title('log 算子检测');
subplot(2,3,5), imshow(BW5);
title('canny 算子检测')
5、实验运行结果:
6、实验问题小结:
在实验验收中,老师问了一个问题,关于中值程序for i=1:p(1)-n+1; for j=1:p(2)-n+1 ;c=x1(i:i+(n-1),j:j+(n-1)); 这段程序的问题,由于当时紧张不知道怎么表达。其实这段程序是为了取出从(i,j)开始的n行n列元素,而i和j的之所以这样取值是因为需要确保确定开始值大小以便能从p范围图像中顺利取出n个元素来。其余问题都已解决,整体较为顺利。谢谢老师耐心的指导。
六、
1、实验内容:
用如下矩阵对图像进行滤波
2、实验目的:
进一步了解图像矩阵滤波的基本知识;
学会用matlab处理图像,实现相应滤波功能。
3.实验代码:
(1)
clear;clc;
I = imread('hat.jpg');
I=I(:,:,1);
I=double(I);
x(:,:,1)=[-1,-1,-1;-1,8,-1;-1,-1,-1];
x(:,:,2)=[-1,-1,-1;2,2,2;-1,-1,-1];
x(:,:,3)=[-1,2,-1;-1,2,-1;-1,2,-1];
x(:,:,4)=[2,-1,-1;-1,2,-1;-1,-1,2];
x(:,:,5)=[-1,-1,2;-1,2,-1;2,-1,-1];
x(:,:,6)=[1,2,1;0,0,0;-1,-2,-1];
x(:,:,7)=[1,0,-1;2,0,-2;1,0,-1];
for k=1:7
for m=2:219
for n=2:219
x1=I(m-1:m+1,n-1:n+1);
y=x(:,:,k).*x1;
y=[y(1,:),y(2,:),y(3,:)];
d(m,n)=sum(y);
end
end
subplot(2,4,k)
imshow(d,[]);title(k);
end
4、实验运行结果:
5、实验小结:
本题的方法有很多,可以用函数调用等方法,大致思路只要正确,最终结果都会正确得出。
七、
实验最后总结:
通过总共两次实验,多题的训练,自己对数字图像有关知识以及matlab的运用有了进一步的了解,每一次实验的问题总结都在每一题后写出。
总的来说,在本课程实验中,平时上机遇到许多问题,当然也收获许多,在验收过程中自己也很紧张,可能语言表达不是太好,这方面依然需要锻炼,做出这些题目不但要掌握课上老师所讲的知识点,也要查阅许多资料,我自己查了很多资料,参考许多东西,然后经过自己的思考,动手操作改进编写,然后变成自己的东西,学会其中的运作思想,我想这也就是学习的很重要的吧,当然,平时也要细心,可能问题就是一点,需要耐心的去排除错误,也要不能畏难,遇到困难也要努力解决,吃透知识然后向老师请教,验收。不能因为老师说不验收也给过就产生松懈,坚持做下去,一定会收获很多,我想这是我体会最深的一点,老师有耐心,很负责,,也让自己更好的学习理解运用所学的知识,当然自己还有许多做的不够好的地方,我会继续努力学习解决问题。最后再次感谢老师的悉心指导。