数字图像处理理工.doc

数字图像处理初步实验
第10章 数字图像处理初步实验
10.1 数字图像直方图统计
10.1.1实验目的
n 掌握直方图统计的原理和程序设计方法
n 了解各种图像的直方图统计意义及其在实际中的运用。
10.1.2 实验设备
n PC兼容机一台；操作系统为Windows2000 (或WindowsNT、Windows98、WindowsXP）
n 计算机安装CCS 5000或CCS 3.1。
10.1.3实验内容
输入一幅图像，使用MATLAB统计其直方图。然后编写DSP程序，输入同一幅图像,使用所编写的程序统计直方图，比较两个直方图的差异。
10.1.4实验原理
灰度直方图是数字图像处理中一个最简单、最有用的工具，它描述了一幅图像的灰度级内容。任何一幅图像的直方图都包括了可观的信息，给出了一幅图像概貌的描述，某些类型的图像可由其直方图完全描述。在数字图像ƒ(m,n)情况下，灰度直方图是灰度值的函数，描述的是图像中具有该灰度值的像素个数，其横坐标表示像素的灰度级别，纵坐标是该灰度出现的频率(像素个数与图像像素总数之比)。设图像像素的灰度值为r0，r1，…，rL-1，则概率Ρ(ri)为：
（i=0，1，…，L-1） （10.1）
且有

尽管灰度直方图不能表示出某灰度值在什么位置，更不能直接显示图像内容，但是具有统计特性的直方图却能描述该图像的灰度分布特性，使人们从中得到诸如总体明亮程度、对比度和对象物的可分性等与图像质量有关的灰度分布概况，成为一些处理方法的重要依据。
在MATLAB中，imhist函数可以用来计算和绘制图像的直方图。对于lena图像，要想显示其直方图可以用下面的程序来完成，其结果见图10-2。
I = imread('lena.bmp'); //输入图像
imshow(I)
figure,imhist(I,64) //直方图统计

图10-1 原始图像及其灰度直方图
10.1.5实验步骤
1．实验准备
设置软件仿真模式，启动CCS。
2．建立工程
建立两个C语言源文件、建立一个命令文件，然后把这三个文件添加到工程。在每个文件中输入参考程序中的相应内容，就得到工程界面如图10-2所示。

图10-2 工程界面
3．编译并装载程序
4．设置原始图像观察窗口
选择菜单View→Graph→Image，作如下设置如图10-3所示。

图10-3 设置原始图像观察窗口
5．设置直方图观察窗口
选择菜单“View→Graph→Time/Frequency”，作如下设置如图10-4所示。

图10-4 设置直方图观察窗口
6．设置断点
在程序Hist.c中，有注释断点的语句上设置软件断点。
7．运行程序
全速运行到各个断点，观察图像和直方图统计结果，如图10-5所示。

图10-5 原始图像及直方图统计所得的图像
10.1.6程序清单参考
1. 主程序Hist.c
#define IMAGEWIDTH 80 //图像的宽度
#define IMAGEHEIGHT 80 //图像的高度
#define MODEPHOTO1 1
#define MODEPHOTO2 2
void InitImage(unsigned int nMode,unsigned char *pImage,int nWidth,int nHeight);
void Histogram(unsigned char *pImage,int nWidth,int nHeight,float fHisto[256]);
unsigned char dbImage[IMAGEWIDTH*IMAGEHEIGHT];
float fHistogram[256];
int main()
{
InitImage(MODEPHOTO1,dbImage,IMAGEWIDTH,IMAGEHEIGHT); //断点
Histogram(dbImage,IMAGEWIDTH,IMAGEHEIGHT,fHistogram);
InitImage(MODEPHOTO2,dbImage,IMAGEWIDTH,IMAGEHEIGHT); //断点
Histogram(dbImage,IMAGEWIDTH,IMAGEHEIGHT,fHistogram);
while (1); //断点
}
/*以下是直方图统计子程序*/
void Histogram(unsigned char *pImage,int nWidth,int nHeight,float fHisto[256])
{
int i,j;
unsigned int uWork;
unsigned char *pWork;
for ( i=0;i<256;i++ ) fHisto[i]=0.0f;
pWork=pImage;
for ( i=0;i00;
}
}

2. 初始化和读图像子程序Initimage.c
#include
#define MODEPHOTO1 1
#define MODEPHOTO2 2
#define GRAYBARLEVEL 16
void ReadImage(unsigned char *pImage,char *cFileName,int nWidth,int nHeight);
void InitImage(unsigned int nMode,unsigned char *pImage,int nWidth,int nHeight)
{
switch ( nMode )
{
case MODEPHOTO1:
ReadImage(pImage,"..\\lina.bmp",nWidth,nHeight);
break;
case MODEPHOTO2:
ReadImage(pImage,"..\\Couple.bmp",nWidth,nHeight);
break;
default:
break;
}
}
void ReadImage(unsigned char *pImage,char *cFileName,int nWidth,int nHeight)
{
int j;
unsigned char *pWork;
FILE *fp;
if ( fp=fopen(cFileName,"rb" ) )
{
fseek(fp,1078L,SEEK_SET);
pWork=pImage+(nHeight-1)*nWidth;
for ( j=0;jY
{
IRAM: o = 000000c0h l = 00040000h
}
SECTIONS
{
.text > IRAM fill = 0
.stack > IRAM fill = 0
.bss > IRAM fill = 0
.cinit > IRAM fill = 0
.const > IRAM fill = 0
.data > IRAM fill = 0
.far > IRAM fill = 0

.switch > IRAM fill = 0
.sysmem > IRAM fill = 0
.tables > IRAM fill = 0
.cio > IRAM fill = 0
}

10.2 数字图像直方图均衡化增强
10.2.1实验目的
n 掌握直方图均衡化增强的原理和程序设计方法；
n 观察对图像进行直方图均衡化增强的效果。
10.2.2 实验设备
n PC兼容机一台；操作系统为Windows2000 (或WindowsNT、Windows98、WindowsXP）
n 计算机安装CCS 5000或CCS 3.1。
10.2.3实验内容
一幅图像，其动态范围较小而且较暗，反映在直方图上就是其直方图所占据的灰度值范围比较窄，而且集中在低灰度值一边。为了使图像更清晰，我们采用直方图均衡化的方法来增加图像灰度动态范围，增强对比度。
首先根据流程图进行编程，其次输入一幅图像, 将其通过所设计的图像处理系统，再次用CCS的虚拟示波器和频谱仪观察输入/输出图像的变化，以达到对输入图像进行直方图均衡化增强的效果。
10.2.4实验原理
图像增强处理技术一直是图像处理领域中一类非常重要的基本处理技术。通过采取适当的增强处理技术，可以将原本模糊不清甚至根本无法分辨的原始图片，处理成清楚、明晰的富含大量有用信息的可使用图像，因此此类图像处理技术在医学、遥感、微生物、刑侦以及军事等诸多领域得到广泛应用。
直方图均衡化的基本思想是把原始图像的直方图变换成均匀分布的形式，这样就增加了像素灰度值的动态范围，从而达到增强图像整体对比度的效果。具体方法是：
1．给出原始图像的所有灰度级Sk,k=0,1,2…,L-1；
2．统计原始图像各灰度级的像素数 nk；
3．利用直方图计算公式 P(Sk)= nk/n ，（k=0,1,…,L-1），计算原始图像的直方图。Sk 代表图像的第k级灰度值，nk 是图像中具有灰度值Sk的像素个数，n是图像像素的总数；
4.计算原始图像的累计直方图：
（si） 0≤Sk≤1 , k=0,1,2,…,L-1；
5．取整计算： Uk =int[(N-1)tk+k/N]；
6．确定映射关系： ；
7．统计新直方图各灰度级Uk的像素数目nk ；
8．计算新的直方图： pl(tk)=。
其程序流程图如图10-6所示：
开始

结束

用不同参数调用构造图像的函数产生图像
调用直方图统计子程序

调用图像增强子程序
产生增强图像

调用直方图统计
子程序
统计增强图的直方图

直方图统计
子程序
直方图统计结束

初始化直方图统
计数组(赋0 值)

统计图像中各种像素
(取值0-255)的个数
统计直方图

计算各种像素个数占全部图像像素的百分比

图像增强
子程序

图像增强
结束
复制直方图统计数组

按照直方图统计
结果生成均衡化
的灰度增强查找

将原始图像中的像
素逐一按照增强灰
度查找表重新计算
灰度生成增强图

图10-6 直方图均衡化程序流程图
原始图像如图10-7所示，其动态范围较小而且较暗，其直方图如图10-8所示，直方图所占据的灰度值范围比较窄，而且集中在低灰度值一边。为了使图像更清晰，我们采用直方图均衡化的方法来增加图像灰度动态范围，增强对比度。在MATLAB中，可以直接调用J=histeq(I,n)函数来完成这项工作，其中I是原始图像矩阵，J是变换后所得的图像矩阵。用户可以指定均衡后的级数n，默认值为64。
I = imread('tire.tif');
J = histeq(I);
imshow(I)
figure, imshow(J)
figure,imhist(I,64)
figure,imhist(J,64)
处理后的图像和直方图分别如图10-9和图10-10所示。

图10-7 原始图像 图10-8 原始直方图

图10-9 直方图均衡化所得的图像 图10-10 直方图均衡化所得的直方图
从处理前后的图像可以看出，许多在原始图像中看不清的细节，在直方图均衡化处理后所得的图像中都变得十分清晰。
10.2.5实验步骤
1．实验准备
设置软件仿真模式，启动CCS。
2．建立工程
建立两个C语言源文件、建立一个命令文件，然后把这三个文件添加到工程。在每个文件中输入参考程序中的相应内容，就得到工程界面如图10-12所示。

图10-11 工程界面

3．编译并装载程序
4．设置原始图像观察窗口
选择菜单View->Graph->Image，作如下设置如图10-13所示。

图10-12 设置原始图像观察窗口
5．设置处理后图像观察窗口
选择菜单View->Graph->Image，作如下设置如图10-13所示。

图10-13 设置处理后图像观察窗口：
6．设置直方图观察窗口
选择菜单View->Graph->Time/Frequency，作如下设置如图10-14。

图10-14 设置直方图观察窗口
7．设置断点
在程序Hist.c中，有注释断点的语句上设置软件断点。
8．运行程序
全速运行到各个断点，观察图像和直方图均衡化结果，分别如图10-15、图10-16、图10-17和图10-18所示。

图10-15 原始图像及直方图均衡化所得的图像

图10-16 处理前（下半部分）、后（上半部分）的直方图

图10-17 原始图像及直方图均衡化所得的图像

图10-18 处理前（下半部分）、后（上半部分）的直方图

10.2.6程序清单参考
1．主程序hist.c
#define IMAGEWIDTH 80 //图像的宽度
#define IMAGEHEIGHT 80 //图像的高度
#define MODEPHOTO1 1
#define MODEPHOTO2 2
void InitImage(unsigned int nMode,unsigned char *pImage,int nWidth,int nHeight);
void Histogram(unsigned char pImage,int nWidth,int nHeight,float fHisto[256]);
void Enhance(unsigned char *pImage,unsigned char *pImage1,int nWidth,int nHeight,float fHisto[256],float fHisto1[256]);
unsigned char dbImage[IMAGEWIDTH*IMAGEHEIGHT];
unsigned char dbTargetImage[IMAGEWIDTH*IMAGEHEIGHT];
float fHistogram[256],fHistogram1[256];
unsigned char lut[256];
int main()
{
InitImage(MODEPHOTO1,dbImage,IMAGEWIDTH,IMAGEHEIGHT); //断点
Histogram(dbImage,IMAGEWIDTH,IMAGEHEIGHT,fHistogram);
Enhance(dbImage,dbTargetImage,IMAGEWIDTH,IMAGEHEIGHT,fHistogram,fHistogram1);
Histogram(dbTargetImage,IMAGEWIDTH,IMAGEHEIGHT,fHistogram1);
InitImage(MODEPHOTO2,dbImage,IMAGEWIDTH,IMAGEHEIGHT); //断点
Histogram(dbImage,IMAGEWIDTH,IMAGEHEIGHT,fHistogram);
Enhance(dbImage,dbTargetImage,IMAGEWIDTH,IMAGEHEIGHT,fHistogram,fHistogram1);
Histogram(dbTargetImage,IMAGEWIDTH,IMAGEHEIGHT,fHistogram1);
while (1); //断点
}
//直方图统计子程序
void Histogram(unsigned char *pImage,int nWidth,int nHeight,float fHisto[256])
{
int i,j;
unsigned int uWork;
unsigned char *pWork;

for ( i=0;i<256;i++ ) fHisto[i]=0.0f;
pWork=pImage;
for ( i=0;i //图像增强子程序
void Enhance(unsigned char *pImage,unsigned char *pImage1,int nWidth,int nHeight,float fHisto[256],float fHisto1[256])
{
int i,j;
unsigned char *pWork,*pWork1;

for ( i=0;i<256;i++ )
fHisto1[i]=fHisto[i]/100;
for ( i=1;i<256;i++ )
fHisto1[i]+=fHisto1[i-1];
for ( i=0;i<256;i++ )
lut[i]=fHisto1[i]*256;
for ( i=0;i<256;i++ )
if ( lut[i]>=256 )
lut[i]=255;
pWork=pImage; pWork1=pImage1;
for ( i=0;i
#include
#define MODEPHOTO1 1
#define MODEPHOTO2 2
#define GRAYBARLEVEL 16
void ReadImage(unsigned char *pImage,char *cFileName,int nWidth,int nHeight);
void InitImage(unsigned int nMode,unsigned char *pImage,int nWidth,int nHeight)
{
switch ( nMode )
{
case MODEPHOTO1:
ReadImage(pImage,"..\\DSP.bmp",nWidth,nHeight);
break;
case MODEPHOTO2:
ReadImage(pImage,"..\\1.bmp",nWidth,nHeight);
break;
default:
break;
}
}
void ReadImage(unsigned char *pImage,char *cFileName,int nWidth,int nHeight)
{
int j;
unsigned char *pWork;
FILE *fp;
if ( fp=fopen(cFileName,"rb" ) )
{
fseek(fp,1078L,SEEK_SET);
pWork=pImage+(nHeight-1)*nWidth;
for ( j=0;j IRAM: o = 000000c0h l = 00040000h
}
SECTIONS
{
.text > IRAM fill = 0
.stack > IRAM fill = 0
.bss > IRAM fill = 0
.cinit > IRAM fill = 0
.const > IRAM fill = 0
.data > IRAM fill = 0

.far > IRAM fill = 0
.switch > IRAM fill = 0
.sysmem > IRAM fill = 0
.tables > IRAM fill = 0
.cio > IRAM fill = 0
}
10.3数字图像边缘检测(Sobel 算子)
10.3.1实验目的
n 了解边缘检测的算法和用途；
n 学习利用Sobel算子进行边缘检测的程序设计方法；
理解边缘检测在图像分割中的应用；
10.3.2 实验设备
n PC兼容机一台；操作系统为Windows2000 (或WindowsNT、Windows98、WindowsXP）
n 计算机安装CCS 5000或CCS 3.1。
10.3.3实验内容
首先根据流程图进行编程，其次输入一幅图像, 将其通过所设计图像处理系统，再次用CCS的虚拟示波器和频谱仪观察输入/输出图像的变化，以达到对输入图像进行边缘检测的效果。
10.3.4实验原理
边缘检测在图像处理与计算机视觉中占有特殊的位置，它是底层视觉处理中最重要的环节之一。边缘能勾划出目标物体，使观察者一目了然边缘蕴含丰富的内在信息（方向、阶跃性质和形状等）。从本质上说，图像边缘是图像局部特征不连续性（灰度突变、颜色突变和纹理结构突变等）的放映，它标志着一个区域的终结和另一个区域的开始。
边缘提取首先检测出图像局部特征的不连续性，然后再将这些不连续的边缘像素连成完备的边界。边缘的特征是沿边沿走向的像素变化平缓，而垂直于边缘方向的像素变化剧烈。所以，从这个意义上说，提取边缘的算法就是检测出符合边缘特性的边缘像素的算子。目前提取边缘常采用边缘算子法、曲面拟合法和模板匹配法等方法。
两个具有不同灰度值的相邻区域之间总存在边缘，边缘是灰度值不连续的结果，这种不连续性通常可以利用求导数的方法方便地检测到一般常用一阶导数和二阶导数来检测边缘。
边缘检测的基本思想是首先利用边缘增强算子，突出图像中的局部边缘，然后定义像素的“边缘强度”，通过设置门限的方法提取边缘点集。常用的边缘检测算于有Robert算子、Sobel算子、Prewitt以算子、LOG算子和Canny算子，下面以Sobel算子进行介绍。
１．Robert算子
1
0
0
-1
0
-1
1
0
Robert算于是一种利用局部差分算子寻找边缘的算子，其模板如图10-19所示。Robert算子对具有陡峭的低噪声的图像效果较好。

图10-19 Robert算子
在MATLAB 中可以由edge函数实现，语法格式如下：
Bw＝edge(I，‘robert’)
Bw=edge(I，‘robert’，thresh)
[BW，thresh]=edge(I，‘toberts’，…)
BW＝edge(I，'robert')自动选择阈值用Robert算子进行边缘检测。
BW＝edge(I.'robert',thresh)根据所指定的敏感度阀值thresh用Robert算子进行边缘检测，它忽略了所有小于阀值的边缘。当thresh为空（[]）时，自动选择阈值。
[BW，thresh]=edge(I,'robert',…)返回阀值。
２．sobel算子
1
0
-1
2
0
-2
1
0
-1
1
2
1
0
0
0
-1
-2
-1
Sobel算子的两个卷积计算核如图10-20所示，图像中的每个点都用这两个核做卷积，第一个核对垂直边缘响应最大，第二个核对水平边缘响应最大。两个卷积的最大值作为该点的输出值，运算结果是一幅边缘幅度图像。sobel算子对灰度渐变和噪声较多的图像处理得较好。

图10-20 sobel算子
edge函数实现的语法格式如下：
BW＝edge(I，‘'sobel’)
BW=edge(I，‘sobel’，thresh)
Bw=edge(I，‘sobel’thesh,direction)
[BW，thresh]=edge(I，‘sobel’…)
Bw＝edge(I，'sobel')自动选择阀值用Sobel算子进行边缘检测。
BW=edge(I，'sobel',thresh)根据所指定的敏感度阀值thresh，用Sobel算子进行边缘检测，它忽略了所有小于阀值的边缘。当thresh为空时，自动选择阀值。
Bw=edge(I,'sobel',thresh,direction)根据所指定的敏感度阈值thresh，在所指定的方向direction上，用Sobe]算子进行边缘检测。direction可取的字符串值为'horizontal'（水平方向）、'vertical'(垂直方向)或'both'（两个方向）。
[BW,thresh]]=edge(I,'sobel'…)返回阀值。
下面的调用与sobel方法类似。
３．prewitt算子
1
1
1
0
0
0
-1
-1
-1
1
0
-1
1
0
-1
1
0
-1
Prewitt算子的两个卷积计算核如图10-22所示，与使用Sobel算子的方法一样，图像中的每个点都用这两个核做卷积，取最大值作为输出。Prewitt算子也产生一幅边缘幅度图像，也是对灰度渐变和噪声较多的图像处理得较好。

图1021prewitt算子
４．LOG算子
拉普拉斯算子的卷积核如图10-23所示，它是一个二阶算子，将在边缘处产生一个陡峭的零交叉。拉普拉斯算子是一个线性的、移不变的算子，它的传递函数在领域空间的原点是零，因此经拉普拉斯滤波过的图像具有零平均灰度。LOG算子先用高斯低通滤波器将图像进行预先平滑，然后用拉普拉斯算子找出图像中的陡峭边缘，最后用零灰度值进行二值化产生闭合的、连通的轮廓，消除了所有内部点。
-1
-1
-1
-1
8
-1
-1
-1
-1
0
-1
0
-1
4
-1
0
-1
0

图10-22 拉普拉斯算子
５．Canny算子
canny算子检测边缘的方法是寻找图像梯度的局部极大值，梯度是用高斯滤波器的导数计算的。canny方法使用两个阈值来分别检测强边缘和弱边缘，而且仅当弱边缘与强边缘相连时，弱边缘才会包含在输出中。因此，此方法不容易受噪声干扰，能够检测到真正的弱边缘。
６．零交叉方法
零交叉方法先用指定的滤波器对图像进行滤波，然后寻找零交叉点作为边缘。
分别采用上述算子在MATLAB环境下，对图像lena检测边缘的结果如图10-24所示。
I = imread('lena.bmp');
BW1 = edge(I,'sobel');
BW2 = edge(I,'roberts');
BW3 = edge(I,'prewitt');
BW4 = edge(I,'log');
BW5 = edge(I,'canny');
h=fspecial('gaussian',5);
BW6 = edge(I,'zerocross',[],h);
subplot(2,3,1),imshow(BW1)
subplot(2,3,2),imshow(BW2)
subplot(2,3,3),imshow(BW3)
subplot(2,3,4),imshow(BW4)
subplot(2,3,5),imshow(BW5)
subplot(2,3,6),imshow(BW6)

图1023 使用各种边缘检测算子得到的边缘图像

开始

结束
用不同参数调用构造图像的函数产生图像
调用边缘检测子
程序生成检测图

边缘检测
子程序

边缘检测
结束

初始化工作变量

针对图像中每一像素用Sobel 算子进行运算

本实验以Sobel算子为例，程序流程图如图10-24所示：

图10-24 程序流程图
10.3.5实验步骤
1．实验准备
设置软件仿真模式，启动CCS。
2．建立工程
建立两个C语言源文件、建立一个命令文件，然后把这三个文件添加到工程。在每个文件中输入参考程序中的相应内容，就得到工程界面如图10-26所示。

图10-25 工程界面
3．编译并装载程序
4．设置原图像观察窗口
选择菜单View->Graph->Image，作如下设置如图10-26所示。

图10-26 原图像观察窗口设置
5．设置边缘图像观察窗口
选择菜单View->Graph->Image，作如下设置如图10-27所示。

图10-27 边缘图像观察窗口设置
6．设置断点
在程序image.c中，有注释断点的语句上设置软件断点。
7. 运行程序
全速运行到各个断点，观察图像和边缘检测结果，如图10-28、图10-29和图10-30所示。

图10-28 边缘检测前后图像比较

图10-29 边缘检测前后图像比较

图

10-30 边缘检测前后图像比较

10.3.6程序清单参考
1．主程序Image.c
#define IMAGEWIDTH 80 //图像的宽度
#define IMAGEHEIGHT 80 //图像的高度
#define MODEPHOTO1 1
#define MODEPHOTO2 2
void InitImage(unsigned int nMode,unsigned char *pImage,int nWidth,int nHeight);
unsigned char dbImage[IMAGEWIDTH*IMAGEHEIGHT];
unsigned char dbTargetImage[IMAGEWIDTH*IMAGEHEIGHT];
int main()
{
InitImage(MODEPHOTO1,dbImage,IMAGEWIDTH,IMAGEHEIGHT); //断点
Sobel(IMAGEWIDTH,IMAGEHEIGHT);

InitImage(MODEPHOTO2,dbImage,IMAGEWIDTH,IMAGEHEIGHT); //断点
Sobel(IMAGEWIDTH,IMAGEHEIGHT);
while (1); //断点
}
2. 初始化及读取图像Initimage.c
#include
#define MODEPHOTO1 1
#define MODEPHOTO2 2
#define GRAYBARLEVEL 16
void ReadImage(unsigned char *pImage,char *cFileName,int nWidth,int nHeight);
void InitImage(unsigned int nMode,unsigned char *pImage,int nWidth,int nHeight)
{
switch ( nMode )
{
case MODEPHOTO1:
ReadImage(pImage,"..\\DSP.bmp",nWidth,nHeight);
break;
case MODEPHOTO2:
ReadImage(pImage,"..\\1.bmp",nWidth,nHeight);
break;
default:
break;
}
}
void ReadIm