(完整word版)数字图像处理图像编码

数字图像处理上机实习报告（DIP4----DIP7）

学生姓名：杜坤

班级：071123

学号：20121003699

指导老师：傅华明

DIP-4 图像编码

一．题目要求

对图实施费诺-香农编码和解码，计算图像熵，平均码长和冗余度。

二．算法设计

1.测试脚本的程序框图

2.编码程序框图

读入图像的直方图，将图像的灰度值按照概率大小排序，按照香农编码的规则编码。

香农编码将概率由大到小，由上到下排成一排，然后分为两组。是将大的一组概率赋值为0，概率小的一组赋值为1，这是赋值的原则。然后依次的重复，直到每组只有一种输入元素为止。

3.解码程序框图

三．实现代码

1.脚本文件

clear all

load mat

p = impr(a); %统计概率

code = FanoCodeInit(p); %Fano编码初始化

code = FanoEncoder(code);%Fano编码

outstream = FanoCodeStream(a,code); %输出

data = FanoDecoder(outstream,code);%解码

data = reshape(data,8,8); %恢复8*8的形状

data = data'; %转置

I = abs(p.*log2(p));

disp('图像的熵为：');

H = sum(I(:)) %计算熵

disp('图像的平局码长为：')

B = FanoCodeLength(code); %求平均长度

disp('编码冗余度为：');

r = B/H - 1 %求冗余

disp('编码效率为：')

e = H/B %求编码效率

if isequal(a,data)

msgbox('解码后的数据和输入的数据完全吻合');

end

2.统计灰度的概率

function [p]= impr(f)

%概率统计

[m,n] = size(f);

graymax = max(f(:)); %找出灰度最大值，划定统计范围p = zeros(1,graymax + 1);

for i = 1:m

for j = 1:n

x = f(i,j) + 1;

p(x) = p(x) + 1;

end

end

p = p/(m*n);

End

3.码字的初始化

function [code] = FanoCodeInit(p)

%FanoShano码字初始化

[m,n] = size(p);

for i = 1:n

code(i).gray = i - 1;

code(i).p = p(i);

code(i).str = '';

end

%冒泡法排序

for i = 1:n

for j = 1:n-i

if code(j).p > code(j+1).p

temp = code(j);

code(j) = code(j+1);

code(j+1) = temp;

end

end

end

end

4.编码

function [pin] = FanoEncoder(pin)

%FanoShano编码

[m,n] = size(pin);

flag = 1;

while (flag)

start = 1;

stop = 1;

temp = pin(1);

for i = 1:n-1

if isequal(temp.str,pin(i+1).str)

stop = stop + 1;

elseif stop == start

start = i + 1;

stop = start;

temp = pin(i+1);

else

break;

end

end

if stop ~= start

pin = FanoCodeCat(pin,start,stop);

else

if i == n-1

flag = 0; %退出while(flag)的循环

end

end

end

end

5.输出码流

function [outstream] = FanoCodeStream(data,code) [m,n] = size(data);

len = length(code);

outstream = '';

for i = 1:m

for j = 1:n

for k = 1:len

if code(k).gray == data(i,j);

outstream = [outstream,code(k).str];

break;

end

end

end

end

end

6.解码

function [data] = FanoDecoder(instream,code)

len = length(instream);

str = '';gray = 0;flag = 0;data = 0;

for i = 1:len

[gray,flag] = LookUp(code,[str,instream(i)]); if flag

dlen = length(data);

data(dlen+1) = gray;

str = '';

else

str = [str,instream(i)];

end

end

dlen = length(data);

data = data(2:dlen);

end

7.搜索码字

function [data,flag] = LookUp(code,str)

len = length(code);

flag = 0;

data = 0;

for i = 1:len

if isequal(str,code(i).str)

data = code(i).gray;

flag = 1;

break;

end

end

end

8.获得平均码长

function [len_ave] = FanoCodeLength(code)

len = length(code);

len_ave = 0;

for i = 1:len

len_ave = len_ave + code(i).p*length(code(i).str);

end

end

四．结果分析

经过检验之后可以看出，将图像数据进行编码，然后再解码得到的数据和原图像数据完全一致，说明此程序成功编码解码，达到了题目的要求。

算法改进：在编码的时候可以直接将灰度值作为码字的下标，提高编码的效率。

DIP-5 图像分割

一．题目要求

对下图施加高斯噪声，采用LoG 算子对含噪声的图象实施边缘分割，找出该图象的最佳边缘。

二．算法设计

在对图像处理的研究和应用中，人们往往仅对图像中的某些部分感兴趣，这些感兴趣的部分常称为目标或对象，它们一般对应图像中特定的、具有独特性质的区域。图像分割是指根据灰度、彩色、空间纹理、几何形状等特征把图像划分成若干个互不相交的区域，使得这些特征在同一区域内表现出一致性或相似性，而在不同区域间表现出明显的不同，即在一幅图像中把目标从背景中分离出来，以便于进一步处理。图像分割就是指把图像分成互不重叠的区域并提取出感兴趣目标的技术。像的分割有很多种类，边缘分割也有很多种类，LoG 算子是其中一类。

由于在成像时，一个给定像素所对应的场景点，它的周围点对该点的贡献的光强大小呈正态分布，所以平滑函数应能反映不同远近的周围点对给定像素具有不同的平滑作用，因此，平滑函数采用正态分布的高斯函数，即

式中，s 是方差。 用h(x ，y)对图像f(x ，y)的平滑可表示为

g(x ，y)＝h(x ，y)*f(x ，y)

如果令r 是离原点的径向距离，即r2＝x2＋y2，转换，然后对图像g(x ，y)采用拉普拉斯算子进行边缘检测，可得

)

,(*),(),(*e )],(*),([),(2242

2

222

2y x f y x h y x f r y x f y x h y x g r ?=???

?

??--=?=?-

σσσ

2

222e

),(σ

y x y x h +-

=

上式中的h2

?称为高斯—拉普拉斯滤波(Laplacian of Gaussian，LoG)算子，也称为“墨西哥草帽”。它是一个轴对称函数，各向同性，它的一个轴截面如图所示。

由图可见，这个函数在r＝±s处有过零点，在|r|s时为负；可以证明这个算子定义域内的平均值为零，因此将它与图像卷积并不会改变图像的整体动态范围。但由于它相当光滑，因此将它与图像卷积会模糊图像，并且其模糊程度是正比于s的。

?的平滑性质能减少噪声的影响，所以当边缘模糊或噪声较大时，正因为h2

?检测过零点能提供较可靠的边缘位置。在该算子中，s的选择很重要，利用h2

s小时边缘位置精度高，但边缘细节变化多；s大时平滑作用大，但细节损失大，边缘点定位精度低。应根据噪声水平和边缘点定位精度要求适当选取s。 LoG 算子用到的卷积模板一般较大，不过这些模板可以分解为一维卷积来快速计算。通过判断零交叉点及其两侧像素符号的变化来确定边缘点。边缘点两侧的二阶微分是异号的，且正号对应边像点的暗侧，负号对应边像点的亮侧，两侧的符号指示着边缘的起伏走向。

三．实现代码

1.主函数

clc;

i = imread('D:\matlab2011\work\p5-03.tif');

subplot(121);imshow(i);title('原图像');

b=log_edge(i);

subplot(122);imshow(b);title('原图像');

2.LoG算子提取边缘点函数

%下面的代码可以实现LoG算子提取边缘点的功能

function e=log_edge(a)

%该函数实现LoG算子提取边缘点

%输入为图像a，输出为边缘图像e

[m,n]=size(a);

e=repmat(logical(uint8(0)),m,n);

sigma=2;

%产生同样大小的边缘图像e，初始化为0

rr=2:m-1;cc=2:n-1;

fsize=ceil(sigma*3)*2+1;

%选择点数为奇数的滤波器的尺寸fsize>6*sigma；

op=fspecial('log',fsize,sigma);

%产生LoG滤波器

op=op-sum(op(:))/prod(size(op));

%将LoG滤波器的均值变为0

b=filter2(op,a);

%利用LoG算子对图像滤波

thresh=.75*mean2(abs(b(rr,cc)));

%设置过零检测的门限

%寻找滤波后的过零点：+ -和- +表示水平方向从左到右和从右到左过零

%[+ -]'和[- +]'表示垂直方向从上到下和从下到上过零

%这里我们选择边缘点为值为负的点

[rx,cx]=find(b(rr,cc)<0&b(rr,cc+1)>0

&abs(b(rr,cc)-b(rr,cc+1))>thresh); %[- +]的情况

e((rx+1)+cx*m)=1;

[rx,cx]=find(b(rr,cc-1)>0&b(rr,cc)<0

&abs(b(rr,cc-1)-b(rr,cc))>thresh); %[+ -]的情况

e((rx+1)+cx*m)=1;

[rx,cx]=find(b(rr,cc)<0&b(rr+1,cc)>0&abs(b(rr,cc)-b(rr+1,cc))>thresh) ; %[- +]'的情况

e((rx+1)+cx*m)=1;

[rx,cx]=find(b(rr-1,cc)>0&b(rr,cc)<0

&abs(b(rr-1,cc)-b(rr,cc))>thresh); %[+ -]'的情况

e((rx+1)+cx*m)=1;

%某些情况下LoG滤波结果可能正好为0，下面考虑这种情况：

[rz,cz]=find(b(rr,cc)==0);

if~isempty(rz)

%寻找滤波后的过零

%+0-和-0+表示水平方向从左到右和从右到左过零

%[+0-]'和[-0+]'表示垂直方向从上到下和从下到上过零

%边缘正好位于滤波值为零点上

zero=(rz+1)+cz*m; %零点的线性坐标

zz=find(b(zero-1)<0&b(zero+1)>0

&abs(b(zero-1)-b(zero+1))>2*thresh); %[-0+]'情况

e(zero(zz))=1;

zz=find(b(zero-1)>0&b(zero+1)<0&abs(b(zero-1)-b(zero+1))>2*thresh); %[+0-]'情况

e(zero(zz))=1;

zz=find(b(zero-m)<0&b(zero+m)>0

&abs(b(zero-m)-b(zero+m))>2*thresh); %[-0+]情况

e(zero(zz))=1;

zz=find(b(zero-m)>0&b(zero+m)<0

&abs(b(zero-m)-b(zero+m))>2*thresh); %[+0-]情况

e(zero(zz))=1;

end

四．结果分析

用LoG算子进行边缘检测的结果如图。Laplace算子对通过图像进行操作实现边缘检测的时，对离散点和噪声比较敏感。于是，首先对图像进行高斯暖卷积滤波进行降噪处理，再采用Laplace算子进行边缘检测，由上图可见有较好提取边缘的效果。

DIP-6 图像描述

一．题目要求

对下图象提取多边形对象的边缘，根据傅氏描述算子重建该图多边形的边界

二．算法设计

傅立叶描述子（Fourier Descriptor，简称FD）常用来表示单封闭曲线的形状特征，其基本思想是将目标轮廓曲线建模成一维序列，对该序列进行一维的傅立叶变换，从而获得一系列的傅立叶系数，用这些系数对该目标轮廓进行描述。傅立叶描述子方法有一系列优点，如：计算原理简单，描述清晰，具有由粗及精的特性等。计算原理简单可以使得特征提取更加稳定，因为在计算的过程中，无须设置大量控制参数就可以获得结果，计算的一致性好。傅立叶描述子具有明确的物理或几何意义，它比某些特征描述子（如Hu 不变矩）更具直观性。此外，由于任何一个序列经傅立叶变换后，其能量主要集中于少数几个低频傅立叶系数上，因此采用极少的傅立叶系数就可以描绘该序列特征。同时，随着傅立叶系数的增多，该序列的细节特征得以更好地描述。因此，傅立叶描述子对目标轮廓有非常好的由粗及精的描述能力。一个傅立叶描述子的构建包括两步：首先，定义一种好的表示（representation）方法对轮廓曲线进行描述；然后，采用傅立叶

理论对该曲线进行变换。不同的曲线表示方法有不同的特性，一个好的表示方法 应该使最终获得的傅立叶描述子具有尺度、旋转、平移不变性及起始点的无关性。

傅立叶描述子，是物体形状边界曲线的傅立叶变换系数，是物体边界曲线信号的频域分析结果。它是一种描述不受起始点移动尺寸变化及旋转影响的曲线的方法。傅立叶描述子的基本思想，是把坐标的序列点看作复数：

()()jy(k)s k s k =+

即x 轴作为实轴，y 轴作为虚轴，边界的性质不变。这种表示方法的优点，是将 一个二维问题简化成一个一维问题。对s(k)的傅立叶变换为：

1

2/0()()N j uk N k a u s k e π--==∑

傅立叶描述子序列{()}a u 反映了原曲线的形状特征，同时，由于傅立叶变换具有能 量集中性，因此，少量的傅立叶描述子就可以重构出原曲线。

三．实现代码

1.主函数

%图像傅里叶算子边界描述 %生成图像边界矩阵 clc;

r=zeros(15,15); %构造15*15的0矩阵

r(3,8:13)=1; %按要求构造“R”型图像

r(4,7)=1; r(4,13)=1;

r(5,6)=1; r(5,13)=1;

r(6,3:5)=1; r(6,13)=1;

r(7,3)=1; r(7,13)=1;

r(8,3)=1; r(8,13)=1;

r(9,3)=1; r(9,13)=1;

r(10,3)=1; r(10,13)=1;

r(11,3)=1; r(11,12)=1;

r(12,3)=1; r(12,11)=1;

r(13,3:10)=1;

subplot(121);

imshow(r);

title('构建的图像边界');

%主函数

%图像傅里叶算子边界描述

i=1;

for m=1:15

for n=1:15

if r(m,n)==1 %将边界坐标存入S矩阵 s(i,1)=m;

s(i,2)=n;

i=i+1;

end

end

end

z=frdescp(s);

%图像傅里叶算子边界逆描述

z34=ifrdescp(z,34);

z34=uint8(z34);

I=zeros(15,15);

for k=1:34

I(z34(k,1),z34(k,2))=1;

k=k+1;

end

subplot(122);

imshow(I);

title('34个描绘子重构的图像');

2.Frdescp函数

功能：计算边界的傅里叶描绘子s

%傅里叶边界描述算子生成函数

function z=frdescp(s)

[np,nc]=size(s);

if nc~=2

error('S must be of size (np,2)'); %必须为2列的矩阵

end

if np/2~=round(np/2); %若点数不是偶数，则补一个点

s(end+1,:)=s(end,:);

np=np+1;

end

x=0:(np-1);

m=((-1).^x)';

s(:,1)=m.*s(:,1);

s(:,2)=m.*s(:,2);

s=s(:,1)+i*s(:,2);

z=fft(s);

3.Ifrdescp函数

功能：给定一组傅里叶描绘子，该函数可用给定数量的描绘子计算其逆变换，以产生一条封闭的空间曲线

%傅里叶边界描述逆算子函数

function s=ifrdescp(z,nd)

np=length(z);

if nargin==1|| nd>np;

nd=np;

end

x=0:(np-1);

m=((-1).^x)';

d=round((np-nd)/2);

z(1:d)=0;

z(np-d+1:np)=0;

zz=ifft(z);

s(:,1)=real(zz);

s(:,2)=imag(zz);

s(:,1)=m.*s(:,1);

s(:,2)=m.*s(:,2);

四．结果分析

1.34个描绘子重构的图像

2.32个描绘子重构的图像

3.29个描绘子重构的图像

4.25个描绘子重构的图像

以上为用34，,3,29,25个描绘子重建的图像。使用34个描绘子重建的图像显示的边界与原图像相同，随着所用描述子减少，出现失真。本图像较简单，若是复杂图像，随着描述子减少重建图像的变化趋势为：1.显示稍微平滑一些的边界，但产生的形状与原图像十分接近。2.仅保留了边界的主要特征。3.丢失边界的主要特征，出现失真。

描绘算子应该尽可能的对平移、旋转和缩放等改变不敏感。当结果取决于所处理的点的顺序时，要给它们加一个额外的约束，以便使描绘子对起始点不敏感。傅里叶描绘子虽然对几何变化简介不敏感，但这些参数的变化却与描绘子的简单变化有关。

DIP-7 图像匹配

一．题目要求

选择摸板8，将选定模板与10个图实施匹配运算。如果匹配成功，请说明图号和摸板左上角像素在匹配图象所在的坐标。

二．算法分析

设检测对象的模板为t(x，y)，令其中心与图像f(x，y)中的一个像素(i，j)重合，检测t(x，y)和图像重合部分之间的相似度，对图像中所有的像素都进行这样的操作，根据相似度为最大或者超过某一阈值来确定对象物是否

存在。本题目中则用到的是普通的模式匹配，即对每一幅图像的每一个像素进行一一的核对。

三．实现代码

clear all

g = imread('D:\matlab2011\work\m08.tif');

[gm,gn] = size(g);

f(1).num = imread('D:\matlab2011\work\p7-01.tif');

f(2).num = imread('D:\matlab2011\work\p7-02.tif');

f(3).num = imread('D:\matlab2011\work\p7-03.bmp');

f(4).num = imread('D:\matlab2011\work\p7-04.tif');

f(5).num = imread('D:\matlab2011\work\p7-05.tif');

f(6).num = imread('D:\matlab2011\work\p7-06.tif');

f(7).num = imread('D:\matlab2011\work\p7-07.tif');

f(8).num = imread('D:\matlab2011\work\p7-08.tif');

f(9).num = imread('D:\matlab2011\work\p7-09.tif');

f(10).num = imread('D:\matlab2011\work\p7-10.tif');

for len = 1:10

[fm,fn] = size(f(len).num);

for i = 1:fm-gm

for j = 1:fn-gn

flag = 1;

for k = 1:gm

for t = 1:gn

if g(k,t) ~= f(len).num(i+k,j+t);

flag = 0;

break;

end

end

if ~flag

break;

end

end

if flag

break;

end

end

if flag

break;