图像处理 线段 角度

2009-07-23 22:30:00 byxdaz 阅读数 23744
  • 直线类型与线段

    掌握OpenCV核心模块,熟练使用相关API 理解各个API背后的相关算法原理,每个参数意义 有能力解决使用应用场景问题,大量工程代码经验分享 掌握图像处理与视频分析,图像分析与测量编码与开发技巧

    67人学习 贾志刚
    免费试看

图像处理与识别学习小结

 

数字图像处理是对图像进行分析、加工、和处理,使其满足视觉、心理以及其他要求的技术。图像处理是信号处理在图像域上的一个应用。目前大多数的图像是以数字形式存储,因而图像处理很多情况下指数字图像处理。此外,基于光学理论的处理方法依然占有重要的地位。 数字图像处理是信号处理的子类, 另外与计算机科学、人工智能等领域也有密切的关系。 传统的一维信号处理的方法和概念很多仍然可以直接应用在图像处理上,比如降噪、量化等。然而,图像属于二维信号,和一维信号相比,它有自己特殊的一面,处理的方式和角度也有所不同。大多数用于一维信号处理的概念都有其在二维图像信号领域的延伸,它们中的一部分在二维情形下变得十分复杂。同时图像处理也具有自身一些新的概念,例如,连通性、旋转不变性,等等。这些概念仅对二维或更高维的情况下才有非平凡的意义。图像处理中常用到快速傅立叶变换,因为它可以减小数据处理量和处理时间。
数字图像处理应用在以下方面

摄影及印刷 (Photography and printing)

卫星图像处理 (Satellite image processing)

医学图像处理 (Medical image processing)

面孔识别, 特征识别 (Face detection, feature detection, face identification)

显微图像处理 (Microscope image processing)

汽车障碍识别 (Car barrier detection)

 

数字图像基础

图像的基本概念、图像取样和量化、数字图像表示、 空间和灰度级分辨率、图像纹理、像素间的一些基本关系(相邻像素、邻接性、连通性、区域和边界、距离度量)、线性和非线性变换

线性变换:如果变换函数是线性的或是分段线性,这种变换就是线性变换。以线性函数加大图像的对比度的效果是使整幅图像的质量改善。以分段线性函数加大图像中某个(或某几个)亮度区间的对比度的效果是使局部亮度区间的质量得到改善。

非线性变换:当变换函数是非线性时,即为非线性变换。常用的有指数变换和对数变换。

RGB (red green blue): 红绿蓝三基色

CMYK (Cyan-Magenta-Yellow-black inK): 青色-品红-黄色-黑色

HSI (Hue-Saturation-Intensity): 色调-饱和度-强度

DDB (device-dependent bitmap): 设备相关位图

DIB (device-independent bitmap): 设备无关位图

CVBS (Composite Video Broadcast Signal): 复合电视广播信号

YUV(亦称Y Cr Cb)是被欧洲电视系统所采用的一种颜色编码方法(属于PAL制)。

 



 

数字图像存储与显示

图像格式

在计算机中,有两种类型的图:矢量图(vector graphics)和位映象图(bitmapped graphics)。矢量图是用数学方法描述的一系列点、线、弧和其他几何形状,如图(a)所示。因此存放这种图使用的格式称为矢量图格式,存储的数据主要是绘制图形的数学描述;位映象图(bitmapped graphics)也称光栅图(raster graphics),这种图就像电视图像一样,由象点组成的,如图(b),因此存放这种图使用的格式称为位映象图格式,经常简称为位图格式,存储的数据是描述像素的数值。

 

矢量图与位映象图

目前包括bmp格式、gif格式、jpeg格式、jpeg2000格式、tiff格式、psd格式、

Png格式、swf格式、svg格式、pcx格式、dxf格式、wmf格式、emf格式、LIC格式、eps格式、TGA格式。

目前比较出名的图像处理库有很多,比如LEADTOOLSOPENCVLEADTOOLS这个是功能非常强大的图像多媒体库,但是这个是收费注册的。OpenCV 是一个跨平台的中、高层 API 构成,目前包括 300 多个 C 函数。它不依赖与其它的外部库,尽管也可以使用某些外部库。OpenCV 对非商业用途和商业用途都是免费(FREE)的。开源的图像库也有不少,比如:

ImageStoneGIMPCxImage等,虽然它们的功能没有LEADTOOLS强大,但是一般的图像处理是可以应付的。

具体的功能介绍参考:http://blog.csdn.net/byxdaz/archive/2009/03/09/3972293.aspx

OpenCV源代码及文档下载:SOURCEFORGE.NET
http://sourceforge.net/projects/opencvlibrary/

 

 

数字图像增强

图像增强的目的在于改善图像的显示质量,以利于信息的提取和识别。从方法上说,则是设法摒弃一些认为不必要或干扰的信息,而将所需要的信息得以突出出来,以利于分析判读或作进一步的处理。以下介绍几种较为简单的遥感数字图像增强处理方法。

A空间域增强处理

空间域是指图像平面所在的二维空间,空间域图像增强是指在图像平面上应用某种数学模型,通过改变图像像元灰度值达到增强效果,这种增强并不改变像元的位置。空域增强包括空域变换增强与空域滤波增强两种。空域变换增强是基于点处理的增强方法、空域滤波增强是基于邻域处理的增强方法。

1)、空域变换增强

常用的空域变换增强方法包括:对比度增强、直方图增强和图像算术运算等。

对比度增强是一种通过改变图像像元的亮度分布态势,扩展灰度分布区间来改变图像像元对比度,从而改善图像质量的图像处理方法。因为亮度值是辐射强度的反映,所以也称为辐射增强。常用的方法有对比度线性变换和非线性变换。其关键是寻找到一个函数,以此函数对图像中每一个像元进行变换,使像元得到统一的重新分配,构成得到反差增强的图像。

直方图增强

直方图均衡化

     直方图均衡化基本做法是将每个灰度区间等概率分布代替了原来的随机分布,即增强后的图象中每一灰度级的像元数目大致相同。直方图均衡化可使得面积最大的地物细节得以增强,而面积小的地物与其灰度接近的地物进行合并,形成综合地物。减少灰度等级换取对比度的增大。

直方图归一化 

     直方图归一化是把原图像的直方图变换为某种指定形态的直方图或某一参考图像的直方图,然后按着已知的指定形态的直方图调整原图像各像元的灰级,最后得到一个直方图匹配的图像。这种方法主要应用在有一幅很好的图像作为标准的情况下,对另一幅不满意的图像用标准图像的直方图进行匹配处理,以改善被处理图像的质量。如在数字镶嵌时,重叠区影像色调由于时相等原因差异往往很大,利用直方图匹配这一方法后可以改善重叠区影像色调过度,如果镶嵌图像时相相差不大,完全可以作到无缝镶嵌。

数字图像的算术运算

两幅或多幅单波段影像,完成空间配准后,通过一系列运算,可以实现图像增强,达到提取某些信息或去掉某些不必要信息的目的。

  

2)、空域滤波增强

空域变换增强是按像元逐点运算的,从整体上改善图像的质量,并不考虑周围像元影响。空间滤波增强则是以重点突出图像上的某些特征为目的的(如突出边缘或纹理等),通过像元与周围相邻像元的关系,采取空间域中的邻域处理方法进行图像增强。邻域法处理用于去噪声、图像平滑、锐化和相关运算。

图像卷积运算是在空间域上对图像作局部检测的运算,以实现平滑和锐化的目的。具体作法是选定一卷积函数,又称为“M×N窗口模板,如3×35×5等。然后从图像左上角开始开一与模板同样大小的活动窗口,图像窗口与模板像元的亮度值对应相乘再相加。将计算结果赋予中心像元作为其灰度值,然后待移动后重新计算,将计算结果赋予另一个中心像元,以此类推直到全幅图像扫描一遍结束生成新的图像。

平滑是指图像中出现某些亮度变化过大的区域,或出现不该有的亮点(噪声)时,采用平滑方法可以减小变化,使亮度平缓或去掉不必要噪声点。它实际上是使图像中高频成分消退,即平滑图像的细节,降低其反差,保存低频成分,在频域中称为低通滤波。具体方法有:均值平滑、中值滤波、锐化。

锐化的作用在于提高边缘灰度值的变化率,使界线更加清晰。它是增强图像中的高频成分,在频域处理中称为高通滤波,也就是使图像细节的反差提高,也称边缘增强。要突出图像的边缘、线状目标或亮度变化率大的部分常采用锐化方法。一般有三种实现方法:

1)梯度法

    梯度反映了相邻像元的亮度变化率,即图像中如果存在边缘,如湖泊、河流的边界,山脉和道路等,则边缘处有较大的梯度值。对于亮度值较平滑的部分,亮度梯度值较小。因此,找到梯度较大的位置,也就找到边缘,然后再用不同的梯度计算值代替边缘处像元的值,也就突出了边缘,实现了图像的锐化。通常有罗伯特梯度和索伯尔梯度方法。

2)拉普拉斯算法

    拉普拉斯算法的意义与梯度法不同,它不检测均匀的亮度变化,而是检测变化率的变化率,相当于二阶微分。计算出的图像更加突出亮度值突变的位置。

3)定向检测

    当有目的地检测某一方向的边、线或纹理特征时,可选择特定的模板卷积运算作定向检测。可以检测垂直边界、水平边界和对角线边界等,各使用的模板不同

 

B频率域图像增强处理
频域增强指在图像的频率域内,对图像的变换系数(频率成分)直接进行运算,然后通过Fourier逆变换以获得图像的增强效果。

一般来说,图像的边缘和噪声对应Fourier变换中的高频部分,所以低通滤波能够平滑图像、去除噪声。

图像灰度发生聚变的部分与频谱的高频分量对应,所以采用高频滤波器衰减或抑制低频分量,能够对图像进行锐化处理。

频域,就是由图像f(x,y)的二维傅立叶变换和相应的频率变量(u,v)的值所组成的空间。在空间域图像强度的变化模式(或规律)可以直接在该空间得到反应。F(0,0)是频域中的原点,反应图像的平均灰度级,即图像中的直流成分;低频反映图像灰度发生缓慢变化的部分;而高频对应图像中灰度发生更快速变化的部分,如边缘、噪声等。但频域不能反应图像的空间信息。

 

 

二维DFT及其反变换、Fast FT

关于这方面的内容需要参考数学知识。

空域和频域滤波间的对应关系:

卷积定理是空域和频域滤波的最基本联系纽带。二维卷积定理:

 

 

 


基本计算过程:

  1. 取函数h(m,n)关于原点的镜像,得到h(-m,-n)
  2. 对某个(x,y),使h(-m,-n)移动相应的距离,得到h(x-m,y-n)
  3. 对积函数f(m,n)h(x-m,y-n)(m,n)的取值范围内求和
  4. 位移是整数增量,对所有的(x,y)重复上面的过程,直到两个函数:f(m,n)h(x-m,y-n)不再有重叠的部分。

 

傅立叶变换是空域和频域的桥梁,关于两个域滤波的傅立叶变换对:

 

 

 

 

 

 

 

 


频域与空域滤波的比较:

1. 对具有同样大小的空域和频率滤波器:h(x,y), H(u,v),频域计算(由于FFT)往往更有效(尤其是图像尺寸比较大时)。但对在空域中用尺寸较小的模板就能解决的问题,则往往在空域中直接操作。

2. 频域滤波虽然更直接,但如果可以使用较小的滤波器,还是在空域计算为好。    因为省去了计算傅立叶变换及反变换等步骤。

3. 由于更多的直观性,频率滤波器设计往往作为空域滤波器设计的向导。

 

平滑的频率域滤波器类型
、理想低通滤波器
、巴特沃思低通滤波器
、高斯低通滤波器
频率域锐化滤波器类型
理想高通滤波器
巴特沃思高通滤波器

高斯型高通滤波器

频率域的拉普拉斯算子
钝化模板、高频提升滤波和高频加强滤波
频率域图像增强处理的过程:

 

 

图像复原
图像复原:试图利用退化过程的先验知识,去恢复已被退化图像的本来面目。

 

图像复原的基本思路:先建立退化的数学模型,然后根据该模型对退化图像进行拟合。

图像复原模型可以用连续数学和离散数学处理,处理项的实现可在空间域卷积,或在频域相乘。 
参考资料:
http://download.csdn.net/source/1513324

 


边缘检测

数字图像的边缘检测是图像分割、目标区域的识别、区域形状提取等图像分析领域十分重要的基础,图像理解和分析的第一步往往就是边缘检测,目前它以成为机器视觉研究领域最活跃的课题之一,在工程应用中占有十分重要的地位。所谓边缘就是指图像局部亮度变化最显著的部分,它是检测图像局部变化显著变化的最基本的运算。边缘的记录有链码表和线段表2种,链码表适合计算周长,线段表容易计算面积以及相关的,他们之间可以相互的转换

常见的边缘检测算法:

Roberts边缘检测算子

Sobel边缘算子

Prewitt边缘算子

Kirsch边缘算子

CANNY边缘检测

 


图像压缩
图像压缩是数据压缩技术在数字图像上的应用,它的目的是减少图像数据中的冗余信息从而用更加高效的格式存储和传输数据。图像压缩可以是有损数据压缩也可以是无损数据压缩。对于如绘制的技术图、图表或者漫画优先使用无损压缩,这是因为有损压缩方法,尤其是在低的位速条件下将会带来压缩失真。如医疗图像或者用于存档的扫描图像等这些有价值的内容的压缩也尽量选择无损压缩方法。有损方法非常适合于自然的图像,例如一些应用中图像的微小损失是可以接受的(有时是无法感知的),这样就可以大幅度地减小位速。

无损图像压缩方法有:

行程长度编码

熵编码法

LZW算法

有损压缩方法有:

将色彩空间化减到图像中常用的颜色。所选择的颜色定义在压缩图像头的调色板中,图像中的每个像素都用调色板中颜色索引表示。这种方法可以与 抖动(en:dithering)一起使用以模糊颜色边界。

色度抽样,这利用了人眼对于亮度变化的敏感性远大于颜色变化,这样就可以将图像中的颜色信息减少一半甚至更多。

变换编码,这是最常用的方法。首先使用如离散余弦变换(DCT)或者小波变换这样的傅立叶相关变换,然后进行量化和用熵编码法压缩。

分形压缩(en:Fractal compression)。



形态学图像处理
 
膨胀与腐蚀

 膨胀
腐蚀
开操作与闭操作
击中或击不中变换
一些基本的形态学算法

边界提取
区域填充
连通分量的提取
凸壳
细化
粗化
骨架

裁剪


图像分割
图像分割是指通过某种方法,使得画面场景中的目标物被分为不同的类别。通常图像分割的实现方法是,将图像分为“黑”、“白”两类,这两类分别代表了两个不同的对象。

图像分割方法:阈值分割区域分割、数学形态学、模式识别方法

A、阈值分割包括以下几种:

1)由直方图灰度分布选择阈值

2)双峰法选择阈值

3)迭代法选取阈值

     原理如下,很好理解。

     迭代法是基于逼近的思想,其步骤如下:
      1. 求出图象的最大灰度值和最小灰度值,分别记为ZMAX和ZMIN,令初始阈值T0=(ZMAX+ZMIN)/2;
     2. 根据阈值TK将图象分割为前景和背景,分别求出两者的平均灰度值ZO和ZB;
     3. 求出新阈值TK+1=(ZO+ZB)/2;
     4. 若TK=TK+1,则所得即为阈值;否则转2,迭代计算。

4 )大津法选择阈值

大津法是属于最大类间方差法,它是自适应计算单阈值的简单高效方法,或者叫(Otsu

大津法由大津于1979年提出,对图像Image,记t为前景与背景的分割阈值,前景点数占图像比例为w0,平均灰度为u0;背景点数占图像比例为w1,平均灰度为u1。图像的总平均灰度为:u=w0*u0+w1*u1。从最小灰度值到最大灰度值遍历t,当t使得值g=w0*(u0-u)2+w1*(u1-u)2 最大时t即为分割的最佳阈值。对大津法可作如下理解:该式实际上就是类间方差值,阈值t分割出的前景和背景两部分构成了整幅图像,而前景取值u0,概率为 w0,背景取值u1,概率为w1,总均值为u,根据方差的定义即得该式。因方差是灰度分布均匀性的一种度量,方差值越大,说明构成图像的两部分差别越大, 当部分目标错分为背景或部分背景错分为目标都会导致两部分差别变小,因此使类间方差最大的分割意味着错分概率最小。直接应用大津法计算量较大,因此一般采用了等价的公式g=w0*w1*(u0-u1)2

5)由灰度拉伸选择阈值

大津法是较通用的方法,但是它对两群物体在灰度不明显的情况下会丢失一些整体信息。因此为了解决这种现象采用灰度拉伸的增强大津法。在大津法的思想上增加灰度的级数来增强前两群物体的灰度差。对于原来的灰度级乘上同一个系数,从而扩大了图像灰度的级数。试验结果表明不同的拉伸系数,分割效果差别比较大。

 

B、区域的分割

区域生长、区域分离与合并

 区域生长算法


C基于形态学分水岭的分割

分水岭分割算法


图像特征提取与匹配

常用的图像特征有颜色特征、纹理特征、形状特征、空间关系特征。

A 颜色特征

特点:颜色特征是一种全局特征,描述了图像或图像区域所对应的景物的表面性质。一般颜色特征是基于像素点的特征,此时所有属于图像或图像区域的像素都有各自的贡献。由于颜色对图像或图像区域的方向、大小等变化不敏感,所以颜色特征不能很好地捕捉图像中对象的局部特征。

常用的特征提取与匹配方法:

颜色直方图

其优点在于:它能简单描述一幅图像中颜色的全局分布,即不同色彩在整幅图像中所占的比例,特别适用于描述那些难以自动分割的图像和不需要考虑物体空间位置的图像。其缺点在于:它无法描述图像中颜色的局部分布及每种色彩所处的空间位置,即无法描述图像中的某一具体的对象或物体。

颜色直方图特征匹配方法:直方图相交法、距离法、中心距法、参考颜色表法、累加颜色直方图法。

 

B 纹理特征

纹理特征的提取方法比较简单,它是用一个活动的窗口在图像上连续滑动,分别计算出窗口中的方差、均值、最大值、最小值及二者之差和信息熵等,

形成相应的纹理图像,当目标的光谱特性比较接近时,纹理特征对于区分目标可以起到积极的作用。选取适当的数据动态变化范围,进行纹理特征提取后,使影像的纹理特征得到突出,有利于提取构造信息。

特点:纹理特征也是一种全局特征,它也描述了图像或图像区域所对应景物的表面性质。但由于纹理只是一种物体表面的特性,并不能完全反映出物体的本质属性,所以仅仅利用纹理特征是无法获得高层次图像内容的。与颜色特征不同,纹理特征不是基于像素点的特征,它需要在包含多个像素点的区域中进行统计计算。在模式匹配中,这种区域性的特征具有较大的优越性,不会由于局部的偏差而无法匹配成功。作为一种统计特征,纹理特征常具有旋转不变性,并且对于噪声有较强的抵抗能力。但是,纹理特征也有其缺点,一个很明显的缺点是当图像的分辨率变化的时候,所计算出来的纹理可能会有较大偏差。另外,由于有可能受到光照、反射情况的影响,从2-D图像中反映出来的纹理不一定是3-D物体表面真实的纹理。

常用的特征提取与匹配方法:

纹理特征描述方法分类

1)统计方法统计方法的典型代表是一种称为灰度共生矩阵的纹理特征分析方法Gotlieb Kreyszig 等人在研究共生矩阵中各种统计特征基础上,通过实验,得出灰度共生矩阵的四个关键特征:能量、惯量、熵和相关性。统计方法中另一种典型方法,则是从图像的自相关函数(即图像的能量谱函数)提取纹理特征,即通过对图像的能量谱函数的计算,提取纹理的粗细度及方向性等特征参数

2)几何法

所谓几何法,是建立在纹理基元(基本的纹理元素)理论基础上的一种纹理特征分析方法。纹理基元理论认为,复杂的纹理可以由若干简单的纹理基元以一定的有规律的形式重复排列构成。在几何法中,比较有影响的算法有两种:Voronio 棋盘格特征法和结构法。

3)模型法

模型法以图像的构造模型为基础,采用模型的参数作为纹理特征。典型的方法是随机场模型法,如马尔可夫(Markov)随机场(MRF)模型法和 Gibbs 随机场模型法

4)信号处理法

纹理特征的提取与匹配主要有:灰度共生矩阵、Tamura 纹理特征、自回归纹理模型、小波变换等。

灰度共生矩阵特征提取与匹配主要依赖于能量、惯量、熵和相关性四个参数。Tamura 纹理特征基于人类对纹理的视觉感知心理学研究,提出6种属性,即

:粗糙度、对比度、方向度、线像度、规整度和粗略度。自回归纹理模型(simultaneous auto-regressive, SAR)是马尔可夫随机场(MRF)模型的一种应用实例。

 

C形状特征

特点:各种基于形状特征的检索方法都可以比较有效地利用图像中感兴趣的目标来进行检索,但它们也有一些共同的问题,

常用的特征提取与匹配方法:

通常情况下,形状特征有两类表示方法,一类是轮廓特征,另一类是区域特征。图像的轮廓特征主要针对物体的外边界,而图像的区域特征则关系到整个形状区域。

几种典型的形状特征描述方法:

1)边界特征法该方法通过对边界特征的描述来获取图像的形状参数。其中Hough 变换检测平行直线方法和边界方向直方图方法是经典方法。Hough 变换是利用图像全局特性而将边缘像素连接起来组成区域封闭边界的一种方法,其基本思想是点—线的对偶性;边界方向直方图法首先微分图像求得图像边缘,然后,做出关于边缘大小和方向的直方图,通常的方法是构造图像灰度梯度方向矩阵。

2)傅里叶形状描述符法

傅里叶形状描述符(Fourier shape descriptors)基本思想是用物体边界的傅里叶变换作为形状描述,利用区域边界的封闭性和周期性,将二维问题转化为一维问题。

由边界点导出三种形状表达,分别是曲率函数、质心距离、复坐标函数。

3)几何参数法

形状的表达和匹配采用更为简单的区域特征描述方法,例如采用有关形状定量测度(如矩、面积、周长等)的形状参数法(shape factor)。在 QBIC 系统中,便是利用圆度、偏心率、主轴方向和代数不变矩等几何参数,进行基于形状特征的图像检索。

 

D空间关系特征

特点:所谓空间关系,是指图像中分割出来的多个目标之间的相互的空间位置或相对方向关系,这些关系也可分为连接/邻接关系、交叠/重叠关系和包含/包容关系等。通常空间位置信息可以分为两类:相对空间位置信息和绝对空间位置信息。前一种关系强调的是目标之间的相对情况,如上下左右关系等,后一种关系强调的是目标之间的距离大小以及方位。显而易见,由绝对空间位置可推出相对空间位置,但表达相对空间位置信息常比较简单。
空间关系特征的使用可加强对图像内容的描述区分能力,但空间关系特征常对图像或目标的旋转、反转、尺度变化等比较敏感。另外,实际应用中,仅仅利用空间信息往往是不够的,不能有效准确地表达场景信息。为了检索,除使用空间关系特征外,还需要其它特征来配合。

常用的特征提取与匹配方法:

提取图像空间关系特征可以有两种方法:一种方法是首先对图像进行自动分割,划分出图像中所包含的对象或颜色区域,然后根据这些区域提取图像特征,并建立索引;另一种方法则简单地将图像均匀地划分为若干规则子块,然后对每个图像子块提取特征,并建立索引。

 

 

模式识别

模式识别是一种从大量信息和数据出发,在专家经验和已有认识的基础上,利用计算机和数学推理的方法对形状、模式、曲线、数字、字符格式和图形自动完成识别的过程。模式识别包括相互关联的两个阶段,即学习阶段和实现阶段,前者是对样本进行特征选择,寻找分类的规律,后者是根据分类规律对未知样本集进行分类和识别。广义的模式识别属计算机科学中智能模拟的研究范畴,内容非常广泛,包括声音和语言识别、文字识别、指纹识别、声纳信号和地震信号分析、照片图片分析、化学模式识别等等。计算机模式识别实现了部分脑力劳动自动化。

模式识别--对表征事物或现象的各种形式的(数值的,文字的和逻辑关系的)信息进行处理和分析,以对事物或现象进行描述、辨认、分类和解释的过程,是信息科学和人工智能的重要组成部分。

模式还可分成抽象的和具体的两种形式。前者如意识、思想、议论等,属于概念识别研究的范畴,是人工智能的另一研究分支。我们所指的模式识别主要是对语音波形、地震波、心电图、脑电图、图片、文字、符号、三位物体和景物以及各种可以用物理的、化学的、生物的传感器对对象进行测量的具体模式进行分类和辨识。

模式识别问题指的是对一系列过程或事件的分类与描述,具有某些相类似的性质的过程或事件就分为一类。模式识别问题一般可以应用以下4种方法进行分析处理。

模版比对:

统计模式识别方法:统计模式识别方法是受数学中的决策理论的启发而产生的一种识别方法,它一般假定被识别的对象或经过特征提取向量是符合一定分布规律的随机变量。其基本思想是将特征提取阶段得到的特征向量定义在一个特征空间中,这个空间包含了所有的特征向量,不同的特征向量,或者说不同类别的对象都对应于空间中的一点。在分类阶段,则利用统计决策的原理对特征空间进行划分,从而达到识别不同特征的对象的目的。统计模式识别中个应用的统计决策分类理论相对比较成熟,研究的重点是特征提取。统计模式识别的基本原理是:有相似性的样本在模式空间中互相接近,并形成集团,即物以类聚。其分析方法是根据模式所测得的特征向量Xi=(xi1,xi2,…,xid)T(i=1,2,…,N),将一个给定的模式归入C个类ω1,ω2,…,ωc中,然后根据模式之间的距离函数来判别分类。其中,T表示转置;N为样本点数;d为样本特征数。

统计模式识别的主要方法有:判别函数法,k近邻分类法,非线性映射法,特征分析法,主因子分析法等。

在统计模式识别中,贝叶斯决策规则从理论上解决了最优分类器的设计问题,但其实施却必须首先解决更困难的概率密度估计问题。BP神经网络直接从观测数据(训练样本)学习,是更简便有效的方法,因而获得了广泛的应用,但它是一种启发式技术,缺乏指定工程实践的坚实理论基础。统计推断理论研究所取得的突破性成果导致现代统计学习理论——VC理论的建立,该理论不仅在严格的数学基础上圆满地回答了人工神经网络中出现的理论问题,而且导出了一种新的学习方法——支撑向量机。

 

人工神经网络模式识别:人工神经网络的研究起源于对生物神经系统的研究。人工神经网络区别于其他识别方法的最大特点是它对待识别的对象不要求有太多的分析与了解,具有一定的智能化处理的特点。

句法结构模式识别:又称结构方法或语言学方法。其基本思想是把一个模式描述为较简单的子模式的组合,子模式又可描述为更简单的子模式的组合,最终得到一个树形的结构描述,在底层的最简单的子模式称为模式基元。在句法方法中选取基元的问题相当于在决策理论方法中选取特征的问题。通常要求所选的基元能对模式提供一个紧凑的反映其结构关系的描述,又要易于用非句法方法加以抽取。显然,基元本身不应该含有重要的结构信息。模式以一组基元和它们的组合关系来描述,称为模式描述语句,这相当于在语言中,句子和短语用词组合,词用字符组合一样。基元组合成模式的规则,由所谓语法来指定。一旦基元被鉴别,识别过程可通过句法分析进行,即分析给定的模式语句是否符合指定的语法,满足某类语法的即被分入该类。

在几种算法中,统计模式识别是最经典的分类识别方法,在图像模式识别中有着非常广泛的应用。

 

 

参考书籍:美国 冈萨雷斯 数字图像处理第二版

2016-07-01 08:44:34 qq_32887673 阅读数 9027
  • 直线类型与线段

    掌握OpenCV核心模块,熟练使用相关API 理解各个API背后的相关算法原理,每个参数意义 有能力解决使用应用场景问题,大量工程代码经验分享 掌握图像处理与视频分析,图像分析与测量编码与开发技巧

    67人学习 贾志刚
    免费试看

图像处理与识别学习小结

 

数字图像处理是对图像进行分析、加工、和处理,使其满足视觉、心理以及其他要求的技术。图像处理是信号处理在图像域上的一个应用。目前大多数的图像是以数字形式存储,因而图像处理很多情况下指数字图像处理。此外,基于光学理论的处理方法依然占有重要的地位。 数字图像处理是信号处理的子类, 另外与计算机科学、人工智能等领域也有密切的关系。 传统的一维信号处理的方法和概念很多仍然可以直接应用在图像处理上,比如降噪、量化等。然而,图像属于二维信号,和一维信号相比,它有自己特殊的一面,处理的方式和角度也有所不同。大多数用于一维信号处理的概念都有其在二维图像信号领域的延伸,它们中的一部分在二维情形下变得十分复杂。同时图像处理也具有自身一些新的概念,例如,连通性、旋转不变性,等等。这些概念仅对二维或更高维的情况下才有非平凡的意义。图像处理中常用到快速傅立叶变换,因为它可以减小数据处理量和处理时间。 
数字图像处理应用在以下方面 :

摄影及印刷 (Photography and printing)

卫星图像处理 (Satellite image processing)

医学图像处理 (Medical image processing)

面孔识别, 特征识别 (Face detection, feature detection, face identification)

显微图像处理 (Microscope image processing)

汽车障碍识别 (Car barrier detection)

 

数字图像基础

图像的基本概念、图像取样和量化、数字图像表示、 空间和灰度级分辨率、图像纹理、像素间的一些基本关系(相邻像素、邻接性、连通性、区域和边界、距离度量)、线性和非线性变换。

线性变换:如果变换函数是线性的或是分段线性,这种变换就是线性变换。以线性函数加大图像的对比度的效果是使整幅图像的质量改善。以分段线性函数加大图像中某个(或某几个)亮度区间的对比度的效果是使局部亮度区间的质量得到改善。

非线性变换:当变换函数是非线性时,即为非线性变换。常用的有指数变换和对数变换。

RGB (red green blue): 红绿蓝三基色

CMYK (Cyan-Magenta-Yellow-black inK): 青色-品红-黄色-黑色

HSI (Hue-Saturation-Intensity): 色调-饱和度-强度

DDB (device-dependent bitmap): 设备相关位图

DIB (device-independent bitmap): 设备无关位图

CVBS (Composite Video Broadcast Signal): 复合电视广播信号

YUV(亦称Y Cr Cb)是被欧洲电视系统所采用的一种颜色编码方法(属于PAL制)。

 



 

数字图像存储与显示

图像格式

在计算机中,有两种类型的图:矢量图(vector graphics)和位映象图(bitmapped graphics)。矢量图是用数学方法描述的一系列点、线、弧和其他几何形状,如图(a)所示。因此存放这种图使用的格式称为矢量图格式,存储的数据主要是绘制图形的数学描述;位映象图(bitmapped graphics)也称光栅图(raster graphics),这种图就像电视图像一样,由象点组成的,如图(b),因此存放这种图使用的格式称为位映象图格式,经常简称为位图格式,存储的数据是描述像素的数值。

 

矢量图与位映象图

目前包括bmp格式、gif格式、jpeg格式、jpeg2000格式、tiff格式、psd格式、

Png格式、swf格式、svg格式、pcx格式、dxf格式、wmf格式、emf格式、LIC格式、eps格式、TGA格式。

目前比较出名的图像处理库有很多,比如LEADTOOLS、OPENCV,LEADTOOLS这个是功能非常强大的图像多媒体库,但是这个是收费注册的。OpenCV 是一个跨平台的中、高层 API 构成,目前包括300 多个 C 函数。它不依赖与其它的外部库,尽管也可以使用某些外部库。OpenCV 对非商业用途和商业用途都是免费(FREE)的。开源的图像库也有不少,比如:

ImageStone、GIMP、CxImage等,虽然它们的功能没有LEADTOOLS强大,但是一般的图像处理是可以应付的。

具体的功能介绍参考:http://blog.csdn.net/byxdaz/archive/2009/03/09/3972293.aspx

OpenCV源代码及文档下载:SOURCEFORGE.NET
http://sourceforge.net/projects/opencvlibrary/

 

 

数字图像增强

图像增强的目的在于改善图像的显示质量,以利于信息的提取和识别。从方法上说,则是设法摒弃一些认为不必要或干扰的信息,而将所需要的信息得以突出出来,以利于分析判读或作进一步的处理。以下介绍几种较为简单的遥感数字图像增强处理方法。

A空间域增强处理

空间域是指图像平面所在的二维空间,空间域图像增强是指在图像平面上应用某种数学模型,通过改变图像像元灰度值达到增强效果,这种增强并不改变像元的位置。空域增强包括空域变换增强与空域滤波增强两种。空域变换增强是基于点处理的增强方法、空域滤波增强是基于邻域处理的增强方法。

1)、空域变换增强

常用的空域变换增强方法包括:对比度增强、直方图增强和图像算术运算等。

对比度增强是一种通过改变图像像元的亮度分布态势,扩展灰度分布区间来改变图像像元对比度,从而改善图像质量的图像处理方法。因为亮度值是辐射强度的反映,所以也称为辐射增强。常用的方法有对比度线性变换和非线性变换。其关键是寻找到一个函数,以此函数对图像中每一个像元进行变换,使像元得到统一的重新分配,构成得到反差增强的图像。

直方图增强

直方图均衡化

     直方图均衡化基本做法是将每个灰度区间等概率分布代替了原来的随机分布,即增强后的图象中每一灰度级的像元数目大致相同。直方图均衡化可使得面积最大的地物细节得以增强,而面积小的地物与其灰度接近的地物进行合并,形成综合地物。减少灰度等级换取对比度的增大。

直方图归一化 

     直方图归一化是把原图像的直方图变换为某种指定形态的直方图或某一参考图像的直方图,然后按着已知的指定形态的直方图调整原图像各像元的灰级,最后得到一个直方图匹配的图像。这种方法主要应用在有一幅很好的图像作为标准的情况下,对另一幅不满意的图像用标准图像的直方图进行匹配处理,以改善被处理图像的质量。如在数字镶嵌时,重叠区影像色调由于时相等原因差异往往很大,利用直方图匹配这一方法后可以改善重叠区影像色调过度,如果镶嵌图像时相相差不大,完全可以作到无缝镶嵌。

数字图像的算术运算

两幅或多幅单波段影像,完成空间配准后,通过一系列运算,可以实现图像增强,达到提取某些信息或去掉某些不必要信息的目的。

  

2)、空域滤波增强

空域变换增强是按像元逐点运算的,从整体上改善图像的质量,并不考虑周围像元影响。空间滤波增强则是以重点突出图像上的某些特征为目的的(如突出边缘或纹理等),通过像元与周围相邻像元的关系,采取空间域中的邻域处理方法进行图像增强。邻域法处理用于去噪声、图像平滑、锐化和相关运算。

图像卷积运算是在空间域上对图像作局部检测的运算,以实现平滑和锐化的目的。具体作法是选定一卷积函数,又称为“M×N窗口”或“模板”,如3×3或5×5等。然后从图像左上角开始开一与模板同样大小的活动窗口,图像窗口与模板像元的亮度值对应相乘再相加。将计算结果赋予中心像元作为其灰度值,然后待移动后重新计算,将计算结果赋予另一个中心像元,以此类推直到全幅图像扫描一遍结束生成新的图像。

平滑是指图像中出现某些亮度变化过大的区域,或出现不该有的亮点(“噪声”)时,采用平滑方法可以减小变化,使亮度平缓或去掉不必要“噪声”点。它实际上是使图像中高频成分消退,即平滑图像的细节,降低其反差,保存低频成分,在频域中称为低通滤波。具体方法有:均值平滑、中值滤波、锐化。

锐化的作用在于提高边缘灰度值的变化率,使界线更加清晰。它是增强图像中的高频成分,在频域处理中称为高通滤波,也就是使图像细节的反差提高,也称边缘增强。要突出图像的边缘、线状目标或亮度变化率大的部分常采用锐化方法。一般有三种实现方法:

(1)梯度法

    梯度反映了相邻像元的亮度变化率,即图像中如果存在边缘,如湖泊、河流的边界,山脉和道路等,则边缘处有较大的梯度值。对于亮度值较平滑的部分,亮度梯度值较小。因此,找到梯度较大的位置,也就找到边缘,然后再用不同的梯度计算值代替边缘处像元的值,也就突出了边缘,实现了图像的锐化。通常有罗伯特梯度和索伯尔梯度方法。

(2)拉普拉斯算法

    拉普拉斯算法的意义与梯度法不同,它不检测均匀的亮度变化,而是检测变化率的变化率,相当于二阶微分。计算出的图像更加突出亮度值突变的位置。

(3)定向检测

    当有目的地检测某一方向的边、线或纹理特征时,可选择特定的模板卷积运算作定向检测。可以检测垂直边界、水平边界和对角线边界等,各使用的模板不同

 

B、频率域图像增强处理
频域增强指在图像的频率域内,对图像的变换系数(频率成分)直接进行运算,然后通过Fourier逆变换以获得图像的增强效果。

一般来说,图像的边缘和噪声对应Fourier变换中的高频部分,所以低通滤波能够平滑图像、去除噪声。

图像灰度发生聚变的部分与频谱的高频分量对应,所以采用高频滤波器衰减或抑制低频分量,能够对图像进行锐化处理。

频域,就是由图像f(x,y)的二维傅立叶变换和相应的频率变量(u,v)的值所组成的空间。在空间域图像强度的变化模式(或规律)可以直接在该空间得到反应。F(0,0)是频域中的原点,反应图像的平均灰度级,即图像中的直流成分;低频反映图像灰度发生缓慢变化的部分;而高频对应图像中灰度发生更快速变化的部分,如边缘、噪声等。但频域不能反应图像的空间信息。

 

 

二维DFT及其反变换、Fast FT

关于这方面的内容需要参考数学知识。

空域和频域滤波间的对应关系:

卷积定理是空域和频域滤波的最基本联系纽带。二维卷积定理:

 

 

 


基本计算过程:

  1. 取函数h(m,n)关于原点的镜像,得到h(-m,-n)
  2. 对某个(x,y),使h(-m,-n)移动相应的距离,得到h(x-m,y-n)
  3. 对积函数f(m,n)h(x-m,y-n)在(m,n)的取值范围内求和
  4. 位移是整数增量,对所有的(x,y)重复上面的过程,直到两个函数:f(m,n)和h(x-m,y-n)不再有重叠的部分。

 

傅立叶变换是空域和频域的桥梁,关于两个域滤波的傅立叶变换对:

 

 

 

 

 

 

 

 


频域与空域滤波的比较:

1. 对具有同样大小的空域和频率滤波器:h(x,y), H(u,v),频域计算(由于FFT)往往更有效(尤其是图像尺寸比较大时)。但对在空域中用尺寸较小的模板就能解决的问题,则往往在空域中直接操作。

2. 频域滤波虽然更直接,但如果可以使用较小的滤波器,还是在空域计算为好。    因为省去了计算傅立叶变换及反变换等步骤。

3. 由于更多的直观性,频率滤波器设计往往作为空域滤波器设计的向导。

 

平滑的频率域滤波器类型
1 、理想低通滤波器
2 、巴特沃思低通滤波器
3 、高斯低通滤波器
频率域锐化滤波器类型
1 理想高通滤波器
2 巴特沃思高通滤波器

3 高斯型高通滤波器

4 频率域的拉普拉斯算子
5 钝化模板、高频提升滤波和高频加强滤波
频率域图像增强处理的过程:

 

 

图像复原
图像复原:试图利用退化过程的先验知识,去恢复已被退化图像的本来面目。

 

图像复原的基本思路:先建立退化的数学模型,然后根据该模型对退化图像进行拟合。

图像复原模型可以用连续数学和离散数学处理,处理项的实现可在空间域卷积,或在频域相乘。  
参考资料:
http://download.csdn.net/source/1513324

 


边缘检测

数字图像的边缘检测是图像分割、目标区域的识别、区域形状提取等图像分析领域十分重要的基础,图像理解和分析的第一步往往就是边缘检测,目前它以成为机器视觉研究领域最活跃的课题之一,在工程应用中占有十分重要的地位。所谓边缘就是指图像局部亮度变化最显著的部分,它是检测图像局部变化显著变化的最基本的运算。边缘的记录有链码表和线段表2种,链码表适合计算周长,线段表容易计算面积以及相关的,他们之间可以相互的转换。

常见的边缘检测算法:

Roberts边缘检测算子

Sobel边缘算子

Prewitt边缘算子

Kirsch边缘算子

CANNY边缘检测

 


图像压缩
图像压缩是数据压缩技术在数字图像上的应用,它的目的是减少图像数据中的冗余信息从而用更加高效的格式存储和传输数据。图像压缩可以是有损数据压缩也可以是无损数据压缩。对于如绘制的技术图、图表或者漫画优先使用无损压缩,这是因为有损压缩方法,尤其是在低的位速条件下将会带来压缩失真。如医疗图像或者用于存档的扫描图像等这些有价值的内容的压缩也尽量选择无损压缩方法。有损方法非常适合于自然的图像,例如一些应用中图像的微小损失是可以接受的(有时是无法感知的),这样就可以大幅度地减小位速。

无损图像压缩方法有:

行程长度编码

熵编码法

LZW算法

有损压缩方法有:

将色彩空间化减到图像中常用的颜色。所选择的颜色定义在压缩图像头的调色板中,图像中的每个像素都用调色板中颜色索引表示。这种方法可以与 抖动(en:dithering)一起使用以模糊颜色边界。

色度抽样,这利用了人眼对于亮度变化的敏感性远大于颜色变化,这样就可以将图像中的颜色信息减少一半甚至更多。

变换编码,这是最常用的方法。首先使用如离散余弦变换(DCT)或者小波变换这样的傅立叶相关变换,然后进行量化和用熵编码法压缩。

分形压缩(en:Fractal compression)。



形态学图像处理
 膨胀与腐蚀

 膨胀
腐蚀
开操作与闭操作
击中或击不中变换
一些基本的形态学算法

边界提取
区域填充
连通分量的提取
凸壳
细化
粗化
骨架

裁剪


图像分割
图像分割是指通过某种方法,使得画面场景中的目标物被分为不同的类别。通常图像分割的实现方法是,将图像分为“黑”、“白”两类,这两类分别代表了两个不同的对象。

图像分割方法:阈值分割、区域分割、数学形态学、模式识别方法

A、阈值分割包括以下几种:

(1)由直方图灰度分布选择阈值

(2)双峰法选择阈值

(3)迭代法选取阈值

     原理如下,很好理解。

     迭代法是基于逼近的思想,其步骤如下:
      1. 求出图象的最大灰度值和最小灰度值,分别记为ZMAX和ZMIN,令初始阈值T0=(ZMAX+ZMIN)/2;
     2. 根据阈值TK将图象分割为前景和背景,分别求出两者的平均灰度值ZO和ZB;
     3. 求出新阈值TK+1=(ZO+ZB)/2;
     4. 若TK=TK+1,则所得即为阈值;否则转2,迭代计算。

(4 )大津法选择阈值

大津法是属于最大类间方差法,它是自适应计算单阈值的简单高效方法,或者叫(Otsu)

大津法由大津于1979年提出,对图像Image,记t为前景与背景的分割阈值,前景点数占图像比例为w0,平均灰度为u0;背景点数占图像比例为w1,平均灰度为u1。图像的总平均灰度为:u=w0*u0+w1*u1。从最小灰度值到最大灰度值遍历t,当t使得值g=w0*(u0-u)2+w1*(u1-u)2 最大时t即为分割的最佳阈值。对大津法可作如下理解:该式实际上就是类间方差值,阈值t分割出的前景和背景两部分构成了整幅图像,而前景取值u0,概率为 w0,背景取值u1,概率为w1,总均值为u,根据方差的定义即得该式。因方差是灰度分布均匀性的一种度量,方差值越大,说明构成图像的两部分差别越大, 当部分目标错分为背景或部分背景错分为目标都会导致两部分差别变小,因此使类间方差最大的分割意味着错分概率最小。直接应用大津法计算量较大,因此一般采用了等价的公式g=w0*w1*(u0-u1)2。

(5)由灰度拉伸选择阈值

大津法是较通用的方法,但是它对两群物体在灰度不明显的情况下会丢失一些整体信息。因此为了解决这种现象采用灰度拉伸的增强大津法。在大津法的思想上增加灰度的级数来增强前两群物体的灰度差。对于原来的灰度级乘上同一个系数,从而扩大了图像灰度的级数。试验结果表明不同的拉伸系数,分割效果差别比较大。

 

B、区域的分割

区域生长、区域分离与合并

 区域生长算法


C、基于形态学分水岭的分割

分水岭分割算法


图像特征提取与匹配

常用的图像特征有颜色特征、纹理特征、形状特征、空间关系特征。

A 颜色特征

特点:颜色特征是一种全局特征,描述了图像或图像区域所对应的景物的表面性质。一般颜色特征是基于像素点的特征,此时所有属于图像或图像区域的像素都有各自的贡献。由于颜色对图像或图像区域的方向、大小等变化不敏感,所以颜色特征不能很好地捕捉图像中对象的局部特征。

常用的特征提取与匹配方法:

颜色直方图

其优点在于:它能简单描述一幅图像中颜色的全局分布,即不同色彩在整幅图像中所占的比例,特别适用于描述那些难以自动分割的图像和不需要考虑物体空间位置的图像。其缺点在于:它无法描述图像中颜色的局部分布及每种色彩所处的空间位置,即无法描述图像中的某一具体的对象或物体。

颜色直方图特征匹配方法:直方图相交法、距离法、中心距法、参考颜色表法、累加颜色直方图法。

 

B 纹理特征

纹理特征的提取方法比较简单,它是用一个活动的窗口在图像上连续滑动,分别计算出窗口中的方差、均值、最大值、最小值及二者之差和信息熵等,

形成相应的纹理图像,当目标的光谱特性比较接近时,纹理特征对于区分目标可以起到积极的作用。选取适当的数据动态变化范围,进行纹理特征提取后,使影像的纹理特征得到突出,有利于提取构造信息。

特点:纹理特征也是一种全局特征,它也描述了图像或图像区域所对应景物的表面性质。但由于纹理只是一种物体表面的特性,并不能完全反映出物体的本质属性,所以仅仅利用纹理特征是无法获得高层次图像内容的。与颜色特征不同,纹理特征不是基于像素点的特征,它需要在包含多个像素点的区域中进行统计计算。在模式匹配中,这种区域性的特征具有较大的优越性,不会由于局部的偏差而无法匹配成功。作为一种统计特征,纹理特征常具有旋转不变性,并且对于噪声有较强的抵抗能力。但是,纹理特征也有其缺点,一个很明显的缺点是当图像的分辨率变化的时候,所计算出来的纹理可能会有较大偏差。另外,由于有可能受到光照、反射情况的影响,从2-D图像中反映出来的纹理不一定是3-D物体表面真实的纹理。

常用的特征提取与匹配方法:

纹理特征描述方法分类

(1)统计方法统计方法的典型代表是一种称为灰度共生矩阵的纹理特征分析方法Gotlieb 和 Kreyszig等人在研究共生矩阵中各种统计特征基础上,通过实验,得出灰度共生矩阵的四个关键特征:能量、惯量、熵和相关性。统计方法中另一种典型方法,则是从图像的自相关函数(即图像的能量谱函数)提取纹理特征,即通过对图像的能量谱函数的计算,提取纹理的粗细度及方向性等特征参数

(2)几何法

所谓几何法,是建立在纹理基元(基本的纹理元素)理论基础上的一种纹理特征分析方法。纹理基元理论认为,复杂的纹理可以由若干简单的纹理基元以一定的有规律的形式重复排列构成。在几何法中,比较有影响的算法有两种:Voronio 棋盘格特征法和结构法。

(3)模型法

模型法以图像的构造模型为基础,采用模型的参数作为纹理特征。典型的方法是随机场模型法,如马尔可夫(Markov)随机场(MRF)模型法和 Gibbs 随机场模型法

(4)信号处理法

纹理特征的提取与匹配主要有:灰度共生矩阵、Tamura 纹理特征、自回归纹理模型、小波变换等。

灰度共生矩阵特征提取与匹配主要依赖于能量、惯量、熵和相关性四个参数。Tamura 纹理特征基于人类对纹理的视觉感知心理学研究,提出6种属性,即

:粗糙度、对比度、方向度、线像度、规整度和粗略度。自回归纹理模型(simultaneous auto-regressive, SAR)是马尔可夫随机场(MRF)模型的一种应用实例。

 

C形状特征

特点:各种基于形状特征的检索方法都可以比较有效地利用图像中感兴趣的目标来进行检索,但它们也有一些共同的问题,

常用的特征提取与匹配方法:

通常情况下,形状特征有两类表示方法,一类是轮廓特征,另一类是区域特征。图像的轮廓特征主要针对物体的外边界,而图像的区域特征则关系到整个形状区域。

几种典型的形状特征描述方法:

(1)边界特征法该方法通过对边界特征的描述来获取图像的形状参数。其中Hough 变换检测平行直线方法和边界方向直方图方法是经典方法。Hough 变换是利用图像全局特性而将边缘像素连接起来组成区域封闭边界的一种方法,其基本思想是点—线的对偶性;边界方向直方图法首先微分图像求得图像边缘,然后,做出关于边缘大小和方向的直方图,通常的方法是构造图像灰度梯度方向矩阵。

(2)傅里叶形状描述符法

傅里叶形状描述符(Fourier shape descriptors)基本思想是用物体边界的傅里叶变换作为形状描述,利用区域边界的封闭性和周期性,将二维问题转化为一维问题。

由边界点导出三种形状表达,分别是曲率函数、质心距离、复坐标函数。

(3)几何参数法

形状的表达和匹配采用更为简单的区域特征描述方法,例如采用有关形状定量测度(如矩、面积、周长等)的形状参数法(shape factor)。在 QBIC 系统中,便是利用圆度、偏心率、主轴方向和代数不变矩等几何参数,进行基于形状特征的图像检索。

 

D空间关系特征

特点:所谓空间关系,是指图像中分割出来的多个目标之间的相互的空间位置或相对方向关系,这些关系也可分为连接/邻接关系、交叠/重叠关系和包含/包容关系等。通常空间位置信息可以分为两类:相对空间位置信息和绝对空间位置信息。前一种关系强调的是目标之间的相对情况,如上下左右关系等,后一种关系强调的是目标之间的距离大小以及方位。显而易见,由绝对空间位置可推出相对空间位置,但表达相对空间位置信息常比较简单。
空间关系特征的使用可加强对图像内容的描述区分能力,但空间关系特征常对图像或目标的旋转、反转、尺度变化等比较敏感。另外,实际应用中,仅仅利用空间信息往往是不够的,不能有效准确地表达场景信息。为了检索,除使用空间关系特征外,还需要其它特征来配合。

常用的特征提取与匹配方法:

提取图像空间关系特征可以有两种方法:一种方法是首先对图像进行自动分割,划分出图像中所包含的对象或颜色区域,然后根据这些区域提取图像特征,并建立索引;另一种方法则简单地将图像均匀地划分为若干规则子块,然后对每个图像子块提取特征,并建立索引。

 

 

模式识别

模式识别是一种从大量信息和数据出发,在专家经验和已有认识的基础上,利用计算机和数学推理的方法对形状、模式、曲线、数字、字符格式和图形自动完成识别的过程。模式识别包括相互关联的两个阶段,即学习阶段和实现阶段,前者是对样本进行特征选择,寻找分类的规律,后者是根据分类规律对未知样本集进行分类和识别。广义的模式识别属计算机科学中智能模拟的研究范畴,内容非常广泛,包括声音和语言识别、文字识别、指纹识别、声纳信号和地震信号分析、照片图片分析、化学模式识别等等。计算机模式识别实现了部分脑力劳动自动化。

模式识别--对表征事物或现象的各种形式的(数值的,文字的和逻辑关系的)信息进行处理和分析,以对事物或现象进行描述、辨认、分类和解释的过程,是信息科学和人工智能的重要组成部分。

模式还可分成抽象的和具体的两种形式。前者如意识、思想、议论等,属于概念识别研究的范畴,是人工智能的另一研究分支。我们所指的模式识别主要是对语音波形、地震波、心电图、脑电图、图片、文字、符号、三位物体和景物以及各种可以用物理的、化学的、生物的传感器对对象进行测量的具体模式进行分类和辨识。

模式识别问题指的是对一系列过程或事件的分类与描述,具有某些相类似的性质的过程或事件就分为一类。模式识别问题一般可以应用以下4种方法进行分析处理。

模版比对:

统计模式识别方法:统计模式识别方法是受数学中的决策理论的启发而产生的一种识别方法,它一般假定被识别的对象或经过特征提取向量是符合一定分布规律的随机变量。其基本思想是将特征提取阶段得到的特征向量定义在一个特征空间中,这个空间包含了所有的特征向量,不同的特征向量,或者说不同类别的对象都对应于空间中的一点。在分类阶段,则利用统计决策的原理对特征空间进行划分,从而达到识别不同特征的对象的目的。统计模式识别中个应用的统计决策分类理论相对比较成熟,研究的重点是特征提取。统计模式识别的基本原理是:有相似性的样本在模式空间中互相接近,并形成“集团”,即“物以类聚”。其分析方法是根据模式所测得的特征向量Xi=(xi1,xi2,…,xid)T(i=1,2,…,N),将一个给定的模式归入C个类ω1,ω2,…,ωc中,然后根据模式之间的距离函数来判别分类。其中,T表示转置;N为样本点数;d为样本特征数。

统计模式识别的主要方法有:判别函数法,k近邻分类法,非线性映射法,特征分析法,主因子分析法等。

在统计模式识别中,贝叶斯决策规则从理论上解决了最优分类器的设计问题,但其实施却必须首先解决更困难的概率密度估计问题。BP神经网络直接从观测数据(训练样本)学习,是更简便有效的方法,因而获得了广泛的应用,但它是一种启发式技术,缺乏指定工程实践的坚实理论基础。统计推断理论研究所取得的突破性成果导致现代统计学习理论——VC理论的建立,该理论不仅在严格的数学基础上圆满地回答了人工神经网络中出现的理论问题,而且导出了一种新的学习方法——支撑向量机。

 

人工神经网络模式识别:人工神经网络的研究起源于对生物神经系统的研究。人工神经网络区别于其他识别方法的最大特点是它对待识别的对象不要求有太多的分析与了解,具有一定的智能化处理的特点。

句法结构模式识别:又称结构方法或语言学方法。其基本思想是把一个模式描述为较简单的子模式的组合,子模式又可描述为更简单的子模式的组合,最终得到一个树形的结构描述,在底层的最简单的子模式称为模式基元。在句法方法中选取基元的问题相当于在决策理论方法中选取特征的问题。通常要求所选的基元能对模式提供一个紧凑的反映其结构关系的描述,又要易于用非句法方法加以抽取。显然,基元本身不应该含有重要的结构信息。模式以一组基元和它们的组合关系来描述,称为模式描述语句,这相当于在语言中,句子和短语用词组合,词用字符组合一样。基元组合成模式的规则,由所谓语法来指定。一旦基元被鉴别,识别过程可通过句法分析进行,即分析给定的模式语句是否符合指定的语法,满足某类语法的即被分入该类。

在几种算法中,统计模式识别是最经典的分类识别方法,在图像模式识别中有着非常广泛的应用。

 

 

参考书籍:美国 冈萨雷斯 数字图像处理第二版

2018-06-28 19:08:01 weixin_42052460 阅读数 20932
  • 直线类型与线段

    掌握OpenCV核心模块,熟练使用相关API 理解各个API背后的相关算法原理,每个参数意义 有能力解决使用应用场景问题,大量工程代码经验分享 掌握图像处理与视频分析,图像分析与测量编码与开发技巧

    67人学习 贾志刚
    免费试看

第 1 章 基本的图像操作和处理

本章讲解操作和处理图像的基础知识,将通过大量示例介绍处理图像所需的 Python 工具包,并介绍用于读取图像、图像转换和缩放、计算导数、画图和保存结果等的基本工具。这些工具的使用将贯穿本书的剩余章节。

1.1 PIL:Python图像处理类库

PIL(Python Imaging Library Python,图像处理类库)提供了通用的图像处理功能,以及大量有用的基本图像操作,比如图像缩放、裁剪、旋转、颜色转换等。PIL 是免费的,可以从 http://www.pythonware.com/products/pil/ 下载。

利用 PIL 中的函数,我们可以从大多数图像格式的文件中读取数据,然后写入最常见的图像格式文件中。PIL 中最重要的模块为 Image。要读取一幅图像,可以使用:

from PIL import Image

pil_im = Image.open('empire.jpg')

上述代码的返回值 pil_im 是一个 PIL 图像对象。

图像的颜色转换可以使用 convert() 方法来实现。要读取一幅图像,并将其转换成灰度图像,只需要加上 convert('L'),如下所示:

pil_im = Image.open('empire.jpg').convert('L')

在 PIL 文档中有一些例子,参见 http://www.pythonware.com/library/pil/handbook/index.htm。这些例子的输出结果如图 1-1 所示。

图 1-1:用 PIL 处理图像的例子

1.1.1 转换图像格式

通过 save() 方法,PIL 可以将图像保存成多种格式的文件。下面的例子从文件名列表(filelist)中读取所有的图像文件,并转换成 JPEG 格式:

from PIL import Image
import os

for infile in filelist:
  outfile = os.path.splitext(infile)[0] + ".jpg"
  if infile != outfile:
    try:
      Image.open(infile).save(outfile)
    except IOError:
      print "cannot convert", infile

PIL 的 open() 函数用于创建 PIL 图像对象,save() 方法用于保存图像到具有指定文件名的文件。除了后缀变为“.jpg”,上述代码的新文件名和原文件名相同。PIL 是个足够智能的类库,可以根据文件扩展名来判定图像的格式。PIL 函数会进行简单的检查,如果文件不是 JPEG 格式,会自动将其转换成 JPEG 格式;如果转换失败,它会在控制台输出一条报告失败的消息。

本书会处理大量图像列表。下面将创建一个包含文件夹中所有图像文件的文件名列表。首先新建一个文件,命名为 imtools.py,来存储一些经常使用的图像操作,然后将下面的函数添加进去:

import os
def get_imlist(path):

""" 返回目录中所有JPG 图像的文件名列表"""

return [os.path.join(path,f) for f in os.listdir(path) if f.endswith('.jpg')]

现在,回到 PIL。

1.1.2 创建缩略图

使用 PIL 可以很方便地创建图像的缩略图。thumbnail() 方法接受一个元组参数(该参数指定生成缩略图的大小),然后将图像转换成符合元组参数指定大小的缩略图。例如,创建最长边为 128 像素的缩略图,可以使用下列命令:

pil_im.thumbnail((128,128))

1.1.3 复制和粘贴图像区域

使用 crop() 方法可以从一幅图像中裁剪指定区域:

box = (100,100,400,400)
region = pil_im.crop(box)

该区域使用四元组来指定。四元组的坐标依次是(左,上,右,下)。PIL 中指定坐标系的左上角坐标为(0,0)。我们可以旋转上面代码中获取的区域,然后使用 paste() 方法将该区域放回去,具体实现如下:

region = region.transpose(Image.ROTATE_180)
pil_im.paste(region,box)

1.1.4 调整尺寸和旋转

要调整一幅图像的尺寸,我们可以调用 resize() 方法。该方法的参数是一个元组,用来指定新图像的大小:

out = pil_im.resize((128,128))

要旋转一幅图像,可以使用逆时针方式表示旋转角度,然后调用 rotate() 方法:

out = pil_im.rotate(45)

上述例子的输出结果如图 1-1 所示。最左端是原始图像,然后是灰度图像、粘贴有旋转后裁剪图像的原始图像,最后是缩略图。

1.2 Matplotlib

我们处理数学运算、绘制图表,或者在图像上绘制点、直线和曲线时,Matplotlib 是个很好的类库,具有比 PIL 更强大的绘图功能。Matplotlib 可以绘制出高质量的图表,就像本书中的许多插图一样。Matplotlib 中的 PyLab 接口包含很多方便用户创建图像的函数。Matplotlib 是开源工具,可以从 http://matplotlib.sourceforge.net/ 免费下载。该链接中包含非常详尽的使用说明和教程。下面的例子展示了本书中需要使用的大部分函数。

1.2.1 绘制图像、点和线

尽管 Matplotlib 可以绘制出较好的条形图、饼状图、散点图等,但是对于大多数计算机视觉应用来说,仅仅需要用到几个绘图命令。最重要的是,我们想用点和线来表示一些事物,比如兴趣点、对应点以及检测出的物体。下面是用几个点和一条线绘制图像的例子:

from PIL import Image
from pylab import *

# 读取图像到数组中
im = array(Image.open('empire.jpg'))

# 绘制图像
imshow(im)

# 一些点
x = [100,100,400,400]
y = [200,500,200,500]

# 使用红色星状标记绘制点
plot(x,y,'r*')

# 绘制连接前两个点的线
plot(x[:2],y[:2])

# 添加标题,显示绘制的图像
title('Plotting: "empire.jpg"')
show()

上面的代码首先绘制出原始图像,然后在 x 和 y 列表中给定点的 x 坐标和 y 坐标上绘制出红色星状标记点,最后在两个列表表示的前两个点之间绘制一条线段(默认为蓝色)。该例子的绘制结果如图 1-2 所示。show() 命令首先打开图形用户界面(GUI),然后新建一个图像窗口。该图形用户界面会循环阻断脚本,然后暂停,直到最后一个图像窗口关闭。在每个脚本里,你只能调用一次 show() 命令,而且通常是在脚本的结尾调用。注意,在 PyLab 库中,我们约定图像的左上角为坐标原点。

图像的坐标轴是一个很有用的调试工具;但是,如果你想绘制出较美观的图像,加上下列命令可以使坐标轴不显示:

axis('off')

上面的命令将绘制出如图 1-2 右边所示的图像。

图 1-2:Matplotlib 绘图示例。带有坐标轴和不带坐标轴的包含点和一条线段的图像

在绘图时,有很多选项可以控制图像的颜色和样式。最有用的一些短命令如表 1-1、表 1-2 和表 1-3 所示。使用方法见下面的例子:

plot(x,y)         # 默认为蓝色实线
plot(x,y,'r*')    # 红色星状标记
plot(x,y,'go-')   # 带有圆圈标记的绿线
plot(x,y,'ks:')   # 带有正方形标记的黑色虚线

表1-1:用PyLab库绘图的基本颜色格式命令

颜色

 

'b'

蓝色

'g'

绿色

'r'

红色

'c'

青色

'm'

品红

'y'

黄色

'k'

黑色

'w'

白色

表1-2:用PyLab库绘图的基本线型格式命令

线型

 

'-'

实线

'--'

虚线

':'

点线

表1-3:用PyLab库绘图的基本绘制标记格式命令

标记

 

'.'

'o'

圆圈

's'

正方形

'*'

星形

'+'

加号

'x'

叉号

1.2.2 图像轮廓和直方图

下面来看两个特别的绘图示例:图像的轮廓和直方图。绘制图像的轮廓(或者其他二维函数的等轮廓线)在工作中非常有用。因为绘制轮廓需要对每个坐标 [x, y] 的像素值施加同一个阈值,所以首先需要将图像灰度化:

from PIL import Image
from pylab import *

# 读取图像到数组中
im = array(Image.open('empire.jpg').convert('L'))

# 新建一个图像
figure()
# 不使用颜色信息
gray()
# 在原点的左上角显示轮廓图像
contour(im, origin='image')
axis('equal')
axis('off')

像之前的例子一样,这里用 PIL 的 convert() 方法将图像转换成灰度图像。

图像的直方图用来表征该图像像素值的分布情况。用一定数目的小区间(bin)来指定表征像素值的范围,每个小区间会得到落入该小区间表示范围的像素数目。该(灰度)图像的直方图可以使用 hist() 函数绘制:

figure()
hist(im.flatten(),128)
show()

hist() 函数的第二个参数指定小区间的数目。需要注意的是,因为 hist() 只接受一维数组作为输入,所以我们在绘制图像直方图之前,必须先对图像进行压平处理。flatten() 方法将任意数组按照行优先准则转换成一维数组。图 1-3 为等轮廓线和直方图图像。

图 1-3:用 Matplotlib 绘制图像等轮廓线和直方图

1.2.3 交互式标注

有时用户需要和某些应用交互,例如在一幅图像中标记一些点,或者标注一些训练数据。PyLab 库中的 ginput() 函数就可以实现交互式标注。下面是一个简短的例子:

from PIL import Image
from pylab import *

im = array(Image.open('empire.jpg'))
imshow(im)
print 'Please click 3 points'
x = ginput(3)
print 'you clicked:',x
show()

上面的脚本首先绘制一幅图像,然后等待用户在绘图窗口的图像区域点击三次。程序将这些点击的坐标 [x, y] 自动保存在 x 列表里。

1.3 NumPy

NumPyhttp://www.scipy.org/NumPy/)是非常有名的 Python 科学计算工具包,其中包含了大量有用的思想,比如数组对象(用来表示向量、矩阵、图像等)以及线性代数函数。NumPy 中的数组对象几乎贯穿用于本书的所有例子中 1 数组对象可以帮助你实现数组中重要的操作,比如矩阵乘积、转置、解方程系统、向量乘积和归一化,这为图像变形、对变化进行建模、图像分类、图像聚类等提供了基础。

1PyLab 实际上包含 NumPy 的一些内容,如数组类型。这也是我们能够在 1.2 节使用数组类型的原因。

NumPy 可以从 http://www.scipy.org/Download 免费下载,在线说明文档(http://docs.scipy.org/doc/numpy/)包含了你可能遇到的大多数问题的答案。关于 NumPy 的更多内容,请参考开源书籍 [24]。

1.3.1 图像数组表示

在先前的例子中,当载入图像时,我们通过调用 array() 方法将图像转换成 NumPy 的数组对象,但当时并没有进行详细介绍。NumPy 中的数组对象是多维的,可以用来表示向量、矩阵和图像。一个数组对象很像一个列表(或者是列表的列表),但是数组中所有的元素必须具有相同的数据类型。除非创建数组对象时指定数据类型,否则数据类型会按照数据的类型自动确定。

对于图像数据,下面的例子阐述了这一点:

im = array(Image.open('empire.jpg'))
print im.shape, im.dtype

im = array(Image.open('empire.jpg').convert('L'),'f')
print im.shape, im.dtype

控制台输出结果如下所示:

(800, 569, 3) uint8
(800, 569) float32

每行的第一个元组表示图像数组的大小(行、列、颜色通道),紧接着的字符串表示数组元素的数据类型。因为图像通常被编码成无符号八位整数(uint8),所以在第一种情况下,载入图像并将其转换到数组中,数组的数据类型为“uint8”。在第二种情况下,对图像进行灰度化处理,并且在创建数组时使用额外的参数“f”;该参数将数据类型转换为浮点型。关于更多数据类型选项,可以参考图书 [24]。注意,由于灰度图像没有颜色信息,所以在形状元组中,它只有两个数值。

数组中的元素可以使用下标访问。位于坐标 ij,以及颜色通道 k 的像素值可以像下面这样访问:

value = im[i,j,k]

多个数组元素可以使用数组切片方式访问。切片方式返回的是以指定间隔下标访问该数组的元素值。下面是有关灰度图像的一些例子:

im[i,:] = im[j,:]      # 将第 j 行的数值赋值给第 i 行
im[:,i] = 100          # 将第 i 列的所有数值设为100
im[:100,:50].sum()     # 计算前100 行、前 50 列所有数值的和
im[50:100,50:100]      # 50~100 行,50~100 列(不包括第 100 行和第 100 列)
im[i].mean()           # 第 i 行所有数值的平均值
im[:,-1]               # 最后一列
im[-2,:] (or im[-2])   # 倒数第二行

注意,示例仅仅使用一个下标访问数组。如果仅使用一个下标,则该下标为行下标。注意,在最后几个例子中,负数切片表示从最后一个元素逆向计数。我们将会频繁地使用切片技术访问像素值,这也是一个很重要的思想。

我们有很多操作和方法来处理数组对象。本书将在使用到的地方逐一介绍。你可以查阅在线文档或者开源图书 [24] 获取更多信息。

1.3.2 灰度变换

将图像读入 NumPy 数组对象后,我们可以对它们执行任意数学操作。一个简单的例子就是图像的灰度变换。考虑任意函数 f,它将 0...255 区间(或者 0...1 区间)映射到自身(意思是说,输出区间的范围和输入区间的范围相同)。下面是关于灰度变换的一些例子:

from PIL import Image
from numpy import *

im = array(Image.open('empire.jpg').convert('L'))

im2 = 255 - im # 对图像进行反相处理

im3 = (100.0/255) * im + 100 # 将图像像素值变换到100...200 区间

im4 = 255.0 * (im/255.0)**2 # 对图像像素值求平方后得到的图像

第一个例子将灰度图像进行反相处理;第二个例子将图像的像素值变换到 100...200 区间;第三个例子对图像使用二次函数变换,使较暗的像素值变得更小。图 1-4 为所使用的变换函数图像。图 1-5 是输出的图像结果。你可以使用下面的命令查看图像中的最小和最大像素值:

print int(im.min()), int(im.max())

图 1-4:灰度变换示例。三个例子中所使用函数的图像,其中虚线表示恒等变换

图 1-5:灰度变换。对图像应用图 1-4 中的函数:f(x)=255-x 对图像进行反相处理(左);f(x)=(100/255)x+100 对图像进行变换(中);f(x)=255(x/255)2 对图像做二次变换(右)

如果试着对上面例子查看最小值和最大值,可以得到下面的输出结果:

2 255
0 253
100 200
0 255

array() 变换的相反操作可以使用 PIL 的 fromarray() 函数完成:

pil_im = Image.fromarray(im)

如果你通过一些操作将“uint8”数据类型转换为其他数据类型,比如之前例子中的 im3 或者 im4,那么在创建 PIL 图像之前,需要将数据类型转换回来:

pil_im = Image.fromarray(uint8(im))

如果你并不十分确定输入数据的类型,安全起见,应该先转换回来。注意,NumPy 总是将数组数据类型转换成能够表示数据的“最低”数据类型。对浮点数做乘积或除法操作会使整数类型的数组变成浮点类型。

1.3.3 图像缩放

NumPy 的数组对象是我们处理图像和数据的主要工具。想要对图像进行缩放处理没有现成简单的方法。我们可以使用之前 PIL 对图像对象转换的操作,写一个简单的用于图像缩放的函数。把下面的函数添加到 imtool.py 文件里:

def imresize(im,sz):
  """ 使用PIL 对象重新定义图像数组的大小"""
  pil_im = Image.fromarray(uint8(im))

  return array(pil_im.resize(sz))

我们将会在接下来的内容中使用这个函数。

1.3.4 直方图均衡化

图像灰度变换中一个非常有用的例子就是直方图均衡化。直方图均衡化是指将一幅图像的灰度直方图变平,使变换后的图像中每个灰度值的分布概率都相同。在对图像做进一步处理之前,直方图均衡化通常是对图像灰度值进行归一化的一个非常好的方法,并且可以增强图像的对比度。

在这种情况下,直方图均衡化的变换函数是图像中像素值的累积分布函数(cumulative distribution function,简写为 cdf,将像素值的范围映射到目标范围的归一化操作)。

下面的函数是直方图均衡化的具体实现。将这个函数添加到 imtool.py 里:

def histeq(im,nbr_bins=256):
  """ 对一幅灰度图像进行直方图均衡化"""

  # 计算图像的直方图
  imhist,bins = histogram(im.flatten(),nbr_bins,normed=True)
  cdf = imhist.cumsum() # cumulative distribution function
  cdf = 255 * cdf / cdf[-1] # 归一化

  # 使用累积分布函数的线性插值,计算新的像素值
  im2 = interp(im.flatten(),bins[:-1],cdf)

  return im2.reshape(im.shape), cdf

该函数有两个输入参数,一个是灰度图像,一个是直方图中使用小区间的数目。函数返回直方图均衡化后的图像,以及用来做像素值映射的累积分布函数。注意,函数中使用到累积分布函数的最后一个元素(下标为 -1),目的是将其归一化到 0...1 范围。你可以像下面这样使用该函数:

from PIL import Image
from numpy import *

im = array(Image.open('AquaTermi_lowcontrast.jpg').convert('L'))
im2,cdf = imtools.histeq(im)

图 1-6 和图 1-7 为上面直方图均衡化例子的结果。上面一行显示的分别是直方图均衡化之前和之后的灰度直方图,以及累积概率分布函数映射图像。可以看到,直方图均衡化后图像的对比度增强了,原先图像灰色区域的细节变得清晰。

图 1-6:直方图均衡化示例。左侧为原始图像和直方图,中间图为灰度变换函数,右侧为直方图均衡化后的图像和相应直方图

图 1-7:直方图均衡化示例。左侧为原始图像和直方图,中间图为灰度变换函数,右侧为直方图均衡化后的图像和相应直方图

1.3.5 图像平均

图像平均操作是减少图像噪声的一种简单方式,通常用于艺术特效。我们可以简单地从图像列表中计算出一幅平均图像。假设所有的图像具有相同的大小,我们可以将这些图像简单地相加,然后除以图像的数目,来计算平均图像。下面的函数可以用于计算平均图像,将其添加到 imtool.py 文件里:

def compute_average(imlist):
  """ 计算图像列表的平均图像"""

  # 打开第一幅图像,将其存储在浮点型数组中
  averageim = array(Image.open(imlist[0]), 'f')

  for imname in imlist[1:]:
    try:
      averageim += array(Image.open(imname))
    except:
      print imname + '...skipped'
  averageim /= len(imlist)

  # 返回uint8 类型的平均图像
  return array(averageim, 'uint8')

该函数包括一些基本的异常处理技巧,可以自动跳过不能打开的图像。我们还可以使用 mean() 函数计算平均图像。mean() 函数需要将所有的图像堆积到一个数组中;也就是说,如果有很多图像,该处理方式需要占用很多内存。我们将会在下一节中使用该函数。

1.3.6 图像的主成分分析(PCA)

PCA(Principal Component Analysis,主成分分析)是一个非常有用的降维技巧。它可以在使用尽可能少维数的前提下,尽量多地保持训练数据的信息,在此意义上是一个最佳技巧。即使是一幅 100×100 像素的小灰度图像,也有 10 000 维,可以看成 10 000 维空间中的一个点。一兆像素的图像具有百万维。由于图像具有很高的维数,在许多计算机视觉应用中,我们经常使用降维操作。PCA 产生的投影矩阵可以被视为将原始坐标变换到现有的坐标系,坐标系中的各个坐标按照重要性递减排列。

为了对图像数据进行 PCA 变换,图像需要转换成一维向量表示。我们可以使用 NumPy 类库中的 flatten() 方法进行变换。

将变平的图像堆积起来,我们可以得到一个矩阵,矩阵的一行表示一幅图像。在计算主方向之前,所有的行图像按照平均图像进行了中心化。我们通常使用 SVD(Singular Value Decomposition,奇异值分解)方法来计算主成分;但当矩阵的维数很大时,SVD 的计算非常慢,所以此时通常不使用 SVD 分解。下面就是 PCA 操作的代码:

from PIL import Image
from numpy import *

def pca(X):
  """ 主成分分析:
    输入:矩阵X ,其中该矩阵中存储训练数据,每一行为一条训练数据
    返回:投影矩阵(按照维度的重要性排序)、方差和均值"""

  # 获取维数
  num_data,dim = X.shape

  # 数据中心化
  mean_X = X.mean(axis=0)
  X = X - mean_X

if dim>num_data:
  # PCA- 使用紧致技巧
  M = dot(X,X.T) # 协方差矩阵
  e,EV = linalg.eigh(M) # 特征值和特征向量
  tmp = dot(X.T,EV).T # 这就是紧致技巧
  V = tmp[::-1] # 由于最后的特征向量是我们所需要的,所以需要将其逆转
  S = sqrt(e)[::-1] # 由于特征值是按照递增顺序排列的,所以需要将其逆转
  for i in range(V.shape[1]):
    V[:,i] /= S
else:
  # PCA- 使用SVD 方法
  U,S,V = linalg.svd(X)
  V = V[:num_data] # 仅仅返回前nun_data 维的数据才合理

# 返回投影矩阵、方差和均值
return V,S,mean_X

该函数首先通过减去每一维的均值将数据中心化,然后计算协方差矩阵对应最大特征值的特征向量,此时可以使用简明的技巧或者 SVD 分解。这里我们使用了 range() 函数,该函数的输入参数为一个整数 n,函数返回整数 0...(n-1) 的一个列表。你也可以使用 arange() 函数来返回一个数组,或者使用 xrange() 函数返回一个产生器(可能会提升速度)。我们在本书中贯穿使用 range() 函数。

如果数据个数小于向量的维数,我们不用 SVD 分解,而是计算维数更小的协方差矩阵 XXT 的特征向量。通过仅计算对应前 kk 是降维后的维数)最大特征值的特征向量,可以使上面的 PCA 操作更快。由于篇幅所限,有兴趣的读者可以自行探索。矩阵 V 的每行向量都是正交的,并且包含了训练数据方差依次减少的坐标方向。

我们接下来对字体图像进行 PCA 变换。fontimages.zip 文件包含采用不同字体的字符 a 的缩略图。所有的 2359 种字体可以免费下载 2。假定这些图像的名称保存在列表 imlist 中,跟之前的代码一起保存传在 pca.py 文件中,我们可以使用下面的脚本计算图像的主成分:

2免费字体图像库由 Martin Solli 收集并上传(http://webstaff.itn.liu.se/~marso/)。

from PIL import Image
from numpy import *
from pylab import *
import pca

im = array(Image.open(imlist[0])) # 打开一幅图像,获取其大小
m,n = im.shape[0:2] # 获取图像的大小
imnbr = len(imlist) # 获取图像的数目

# 创建矩阵,保存所有压平后的图像数据
immatrix = array([array(Image.open(im)).flatten()
               for im in imlist],'f')

# 执行 PCA 操作
V,S,immean = pca.pca(immatrix)

# 显示一些图像(均值图像和前 7 个模式)
figure()
gray()
subplot(2,4,1)
imshow(immean.reshape(m,n))
for i in range(7):
  subplot(2,4,i+2)
  imshow(V[i].reshape(m,n))

show()

注意,图像需要从一维表示重新转换成二维图像;可以使用 reshape() 函数。如图 1-8 所示,运行该例子会在一个绘图窗口中显示 8 个图像。这里我们使用了 PyLab 库的 subplot() 函数在一个窗口中放置多个图像。

图 1-8:平均图像(左上)和前 7 个模式(具有最大方差的方向模式)

1.3.7 使用pickle模块

如果想要保存一些结果或者数据以方便后续使用,Python 中的 pickle 模块非常有用。pickle 模块可以接受几乎所有的 Python 对象,并且将其转换成字符串表示,该过程叫做封装(pickling)。从字符串表示中重构该对象,称为拆封(unpickling)。这些字符串表示可以方便地存储和传输。

我们来看一个例子。假设想要保存上一节字体图像的平均图像和主成分,可以这样来完成:

# 保存均值和主成分数据
f = open('font_pca_modes.pkl', 'wb')
pickle.dump(immean,f)
pickle.dump(V,f)
f.close()

在上述例子中,许多对象可以保存到同一个文件中。pickle 模块中有很多不同的协议可以生成 .pkl 文件;如果不确定的话,最好以二进制文件的形式读取和写入。在其他 Python 会话中载入数据,只需要如下使用 load() 方法:

# 载入均值和主成分数据
f = open('font_pca_modes.pkl', 'rb')
immean = pickle.load(f)
V = pickle.load(f)
f.close()

注意,载入对象的顺序必须和先前保存的一样。Python 中有个用 C 语言写的优化版本,叫做 cpickle 模块,该模块和标准 pickle 模块完全兼容。关于 pickle 模块的更多内容,参见 pickle 模块文档页 http://docs.python.org/library/pickle.html

在本书接下来的章节中,我们将使用 with 语句处理文件的读写操作。这是 Python 2.5 引入的思想,可以自动打开和关闭文件(即使在文件打开时发生错误)。下面的例子使用 with() 来实现保存和载入操作:

# 打开文件并保存
with open('font_pca_modes.pkl', 'wb') as f:
  pickle.dump(immean,f)
  pickle.dump(V,f)

# 打开文件并载入
with open('font_pca_modes.pkl', 'rb') as f:
  immean = pickle.load(f)
  V = pickle.load(f)

上面的例子乍看起来可能很奇怪,但 with() 确实是个很有用的思想。如果你不喜欢它,可以使用之前的 open 和 close函数。

作为 pickle 的一种替代方式,NumPy 具有读写文本文件的简单函数。如果数据中不包含复杂的数据结构,比如在一幅图像上点击的点列表,NumPy 的读写函数会很有用。保存一个数组 x 到文件中,可以使用:

savetxt('test.txt',x,'%i')

最后一个参数表示应该使用整数格式。类似地,读取可以使用:

x = loadtxt('test.txt')

你可以从在线文档 http://docs.scipy.org/doc/numpy/reference/generated/numpy.loadtxt.html 了解更多内容。

最后,NumPy 有专门用于保存和载入数组的函数。你可以在上面的在线文档里查看关于 save() 和 load() 的更多内容。

1.4 SciPy

SciPyhttp://scipy.org/) 是建立在 NumPy 基础上,用于数值运算的开源工具包。SciPy 提供很多高效的操作,可以实现数值积分、优化、统计、信号处理,以及对我们来说最重要的图像处理功能。接下来,本节会介绍 SciPy 中大量有用的模块。SciPy 是个开源工具包,可以从 http://scipy.org/Download 下载。

1.4.1 图像模糊

图像的高斯模糊是非常经典的图像卷积例子。本质上,图像模糊就是将(灰度)图像 I 和一个高斯核进行卷积操作:

Iσ = I*

其中 * 表示卷积操作; 是标准差为 σ 的二维高斯核,定义为 :

G_\sigma=\frac{1}{2\pi\sigma}e^{-(x^2+y^2)/2\sigma^2}

高斯模糊通常是其他图像处理操作的一部分,比如图像插值操作、兴趣点计算以及很多其他应用。

SciPy 有用来做滤波操作的 scipy.ndimage.filters 模块。该模块使用快速一维分离的方式来计算卷积。你可以像下面这样来使用它:

from PIL import Image
from numpy import *
from scipy.ndimage import filters

im = array(Image.open('empire.jpg').convert('L'))
im2 = filters.gaussian_filter(im,5)

上面 guassian_filter() 函数的最后一个参数表示标准差。

图 1-9 显示了随着 σ 的增加,一幅图像被模糊的程度。σ 越大,处理后的图像细节丢失越多。如果打算模糊一幅彩色图像,只需简单地对每一个颜色通道进行高斯模糊:

im = array(Image.open('empire.jpg'))
im2 = zeros(im.shape)
for i in range(3):
  im2[:,:,i] = filters.gaussian_filter(im[:,:,i],5)
im2 = uint8(im2)

在上面的脚本中,最后并不总是需要将图像转换成 uint8 格式,这里只是将像素值用八位来表示。我们也可以使用:

im2 = array(im2,'uint8')

来完成转换。

关于该模块更多的内容以及不同参数的选择,请查看 http://docs.scipy.org/doc/scipy/reference/ndimage.html 上 SciPy 文档中的 scipy.ndimage 部分。

图 1-9:使用 scipy.ndimage.filters 模块进行高斯模糊:(a)原始灰度图像;(b)使用 σ=2 的高斯滤波器;(c)使用 σ=5 的高斯滤波器;(d)使用 σ=10 的高斯滤波器

1.4.2 图像导数

整本书中可以看到,在很多应用中图像强度的变化情况是非常重要的信息。强度的变化可以用灰度图像 I(对于彩色图像,通常对每个颜色通道分别计算导数)的 x 和 y 方向导数 Ix 和 Iy 进行描述。

图像的梯度向量为∇I = [IxIy]T。梯度有两个重要的属性,一是梯度的大小

\left|\boldsymbol{\nabla I}\right|=\sqrt{{\boldsymbol{I}_x}^2+{\boldsymbol{I}_y}^2}

它描述了图像强度变化的强弱,一是梯度的角度

α=arctan2(IyIx)

描述了图像中在每个点(像素)上强度变化最大的方向。NumPy 中的 arctan2() 函数返回弧度表示的有符号角度,角度的变化区间为 -π...π。

我们可以用离散近似的方式来计算图像的导数。图像导数大多数可以通过卷积简单地实现:

Ix=I*Dx 和 Iy=I*Dy

对于 Dx 和 Dy,通常选择 Prewitt 滤波器:

D_x=\begin{vmatrix}-1&0&1\\-1&0&1\\ -1&0&1\end{vmatrix} 和 D_y=\begin{vmatrix}-1&-1&-1\\0&0&0\\ 1&1&1\end{vmatrix}

或者 Sobel 滤波器:

D_x=\begin{vmatrix}-1&0&1\\-2&0&2\\ -1&0&1\end{vmatrix} 和 D_y=\begin{vmatrix}-1&-2&-1\\0&0&0\\ 1&2&1\end{vmatrix}

这些导数滤波器可以使用 scipy.ndimage.filters 模块的标准卷积操作来简单地实现,例如:

from PIL import Image
from numpy import *
from scipy.ndimage import filters

im = array(Image.open('empire.jpg').convert('L'))

# Sobel 导数滤波器
imx = zeros(im.shape)
filters.sobel(im,1,imx)

imy = zeros(im.shape)
filters.sobel(im,0,imy)

magnitude = sqrt(imx**2+imy**2)

上面的脚本使用 Sobel 滤波器来计算 x 和 y 的方向导数,以及梯度大小。sobel() 函数的第二个参数表示选择 x 或者 y 方向导数,第三个参数保存输出的变量。图 1-10 显示了用 Sobel 滤波器计算出的导数图像。在两个导数图像中,正导数显示为亮的像素,负导数显示为暗的像素。灰色区域表示导数的值接近于零。

图 1-10:使用 Sobel 导数滤波器计算导数图像:(a)原始灰度图像;(b)x 导数图像;(c)y 导数图像;(d)梯度大小图像

上述计算图像导数的方法有一些缺陷:在该方法中,滤波器的尺度需要随着图像分辨率的变化而变化。为了在图像噪声方面更稳健,以及在任意尺度上计算导数,我们可以使用高斯导数滤波器:

Ix=I*Gσx 和 Iy=I*Gσy

其中,Gσx和 Gσy 表示  在 x 和 y 方向上的导数, 为标准差为 σ 的高斯函数。

我们之前用于模糊的 filters.gaussian_filter() 函数可以接受额外的参数,用来计算高斯导数。可以简单地按照下面的方式来处理:

sigma = 5 # 标准差

imx = zeros(im.shape)
filters.gaussian_filter(im, (sigma,sigma), (0,1), imx)

imy = zeros(im.shape)
filters.gaussian_filter(im, (sigma,sigma), (1,0), imy)

该函数的第三个参数指定对每个方向计算哪种类型的导数,第二个参数为使用的标准差。你可以查看相应文档了解详情。图 1-11 显示了不同尺度下的导数图像和梯度大小。你可以和图 1-9 中做相同尺度模糊的图像做比较。

图 1-11:使用高斯导数计算图像导数:x 导数图像(上),y 导数图像(中),以及梯度大小图像(下);(a)为原始灰度图像,(b)为使用 σ=2 的高斯导数滤波器处理后的图像,(c)为使 用 σ=5 的高斯导数滤波器处理后的图像,(d)为使用 σ=10 的高斯导数滤波器处理后的图像

1.4.3 形态学:对象计数

形态学(或数学形态学)是度量和分析基本形状的图像处理方法的基本框架与集合。形态学通常用于处理二值图像,但是也能够用于灰度图像。二值图像是指图像的每个像素只能取两个值,通常是 0 和 1。二值图像通常是,在计算物体的数目,或者度量其大小时,对一幅图像进行阈值化后的结果。你可以从http://en.wikipedia.org/wiki/Mathematical_morphology 大体了解形态学及其处理图像的方式。

scipy.ndimage 中的 morphology 模块可以实现形态学操作。你可以使用 scipy.ndimage 中的 measurements 模块来实现二值图像的计数和度量功能。下面通过一个简单的例子介绍如何使用它们。

考虑在图 1-12a3 里的二值图像,计算该图像中的对象个数可以通过下面的脚本实现:

3这个图像实际上是图像“分割”后的结果。如果你想知道该图像是如何创建的,可以查看 9.3 节。

from scipy.ndimage import measurements,morphology

# 载入图像,然后使用阈值化操作,以保证处理的图像为二值图像
im = array(Image.open('houses.png').convert('L'))
im = 1*(im<128)

labels, nbr_objects = measurements.label(im)
print "Number of objects:", nbr_objects

上面的脚本首先载入该图像,通过阈值化方式来确保该图像是二值图像。通过和 1 相乘,脚本将布尔数组转换成二进制表示。然后,我们使用 label() 函数寻找单个的物体,并且按照它们属于哪个对象将整数标签给像素赋值。图 1-12b 是labels 数组的图像。图像的灰度值表示对象的标签。可以看到,在一些对象之间有一些小的连接。进行二进制开(binary open)操作,我们可以将其移除:

# 形态学开操作更好地分离各个对象
im_open = morphology.binary_opening(im,ones((9,5)),iterations=2)

labels_open, nbr_objects_open = measurements.label(im_open)
print "Number of objects:", nbr_objects_open

binary_opening() 函数的第二个参数指定一个数组结构元素。该数组表示以一个像素为中心时,使用哪些相邻像素。在这种情况下,我们在 y 方向上使用 9 个像素(上面 4 个像素、像素本身、下面 4 个像素),在 x 方向上使用 5 个像素。你可以指定任意数组为结构元素,数组中的非零元素决定使用哪些相邻像素。参数 iterations 决定执行该操作的次数。你可以尝试使用不同的迭代次数 iterations 值,看一下对象的数目如何变化。你可以在图 1-12c 与图 1-12d 中查看经过开操作后的图像,以及相应的标签图像。正如你想象的一样,binary_closing() 函数实现相反的操作。我们将该函数和在morphology 和 measurements 模块中的其他函数的用法留作练习。你可以从 scipy.ndimage 模块文档http://docs.scipy.org/doc/scipy/reference/ndimage.html 中了解关于这些函数的更多知识。

图 1-12:形态学示例。使用二值开操作将对象分开,然后计算物体的数目:(a)为原始二值图像;(b)为对应原始图像的标签图像,其中灰度值表示物体的标签;(c)为使用开操作后的二值图像;(d)为开操作后图像的标签图像

1.4.4 一些有用的SciPy模块

SciPy 中包含一些用于输入和输出的实用模块。下面介绍其中两个模块:io 和 misc

  1. 读写.mat文件

    如果你有一些数据,或者在网上下载到一些有趣的数据集,这些数据以 Matlab 的 .mat 文件格式存储,那么可以使用 scipy.io 模块进行读取。

    data = scipy.io.loadmat('test.mat')

    上面代码中,data 对象包含一个字典,字典中的键对应于保存在原始 .mat 文件中的变量名。由于这些变量是数组格式的,因此可以很方便地保存到 .mat 文件中。你仅需创建一个字典(其中要包含你想要保存的所有变量),然后使用 savemat() 函数:

    data = {}
    data['x'] = x
    scipy.io.savemat('test.mat',data)

    因为上面的脚本保存的是数组 x,所以当读入到 Matlab 中时,变量的名字仍为 x。关于 scipy.io 模块的更多内容,请参见在线文档 http://docs.scipy.org/doc/scipy/reference/io.html

  2. 以图像形式保存数组

    因为我们需要对图像进行操作,并且需要使用数组对象来做运算,所以将数组直接保存为图像文件 4 非常有用。本书中的很多图像都是这样的创建的。

    imsave() 函数可以从 scipy.misc 模块中载入。要将数组 im 保存到文件中,可以使用下面的命令:

    from scipy.misc import imsave
    imsave('test.jpg',im)

    scipy.misc 模块同样包含了著名的 Lena 测试图像:

    lena = scipy.misc.lena()

    该脚本返回一个 512×512 的灰度图像数组。

4所有 Pylab 图均可保存为多种图像格式,方法是点击图像窗口中的“保存”按钮。

1.5 高级示例:图像去噪

我们通过一个非常实用的例子——图像的去噪——来结束本章。图像去噪是在去除图像噪声的同时,尽可能地保留图像细节和结构的处理技术。我们这里使用 ROF(Rudin-Osher-Fatemi)去噪模型。该模型最早出现在文献 [28] 中。图像去噪对于很多应用来说都非常重要;这些应用范围很广,小到让你的假期照片看起来更漂亮,大到提高卫星图像的质量。ROF 模型具有很好的性质:使处理后的图像更平滑,同时保持图像边缘和结构信息。

ROF 模型的数学基础和处理技巧非常高深,不在本书讲述范围之内。在讲述如何基于 Chambolle 提出的算法 [5] 实现 ROF 求解器之前,本书首先简要介绍一下 ROF 模型。

一幅(灰度)图像 I 的全变差(Total Variation,TV)定义为梯度范数之和。在连续表示的情况下,全变差表示为:

J(\boldsymbol{I})=\int\left|\nabla\boldsymbol{I}\right|\text{dx}            (1.1)

在离散表示的情况下,全变差表示为:

J(\boldsymbol{I})=\sum_{\text{x}}\left|\nabla\boldsymbol{I}\right|

其中,上面的式子是在所有图像坐标 x=[x, y] 上取和。

在 Chambolle 提出的 ROF 模型里,目标函数为寻找降噪后的图像 U,使下式最小:

\min_U\left|\left|\boldsymbol{I}-\boldsymbol{U}\right|\right|^2+2\lambda J(\boldsymbol{U}),

其中范数 ||I-U|| 是去噪后图像 U 和原始图像 I 差异的度量。也就是说,本质上该模型使去噪后的图像像素值“平坦”变化,但是在图像区域的边缘上,允许去噪后的图像像素值“跳跃”变化。

按照论文 [5] 中的算法,我们可以按照下面的代码实现 ROF 模型去噪:

from numpy import *

def denoise(im,U_init,tolerance=0.1,tau=0.125,tv_weight=100):
  """ 使用A. Chambolle(2005)在公式(11)中的计算步骤实现Rudin-Osher-Fatemi(ROF)去噪模型

    输入:含有噪声的输入图像(灰度图像)、U 的初始值、TV 正则项权值、步长、停业条件

    输出:去噪和去除纹理后的图像、纹理残留"""

m,n = im.shape # 噪声图像的大小

# 初始化
U = U_init
Px = im # 对偶域的x 分量
Py = im # 对偶域的y 分量
error = 1

while (error > tolerance):
  Uold = U

  # 原始变量的梯度
  GradUx = roll(U,-1,axis=1)-U # 变量U 梯度的x 分量
  GradUy = roll(U,-1,axis=0)-U # 变量U 梯度的y 分量

  # 更新对偶变量
  PxNew = Px + (tau/tv_weight)*GradUx
  PyNew = Py + (tau/tv_weight)*GradUy
  NormNew = maximum(1,sqrt(PxNew**2+PyNew**2))

  Px = PxNew/NormNew # 更新x 分量(对偶)
  Py = PyNew/NormNew # 更新y 分量(对偶)

  # 更新原始变量
  RxPx = roll(Px,1,axis=1) # 对x 分量进行向右x 轴平移
  RyPy = roll(Py,1,axis=0) # 对y 分量进行向右y 轴平移

  DivP = (Px-RxPx)+(Py-RyPy) # 对偶域的散度
  U = im + tv_weight*DivP # 更新原始变量

  # 更新误差
  error = linalg.norm(U-Uold)/sqrt(n*m);

return U,
2016-06-24 12:55:25 wuxiaobingandbob 阅读数 75694
  • 直线类型与线段

    掌握OpenCV核心模块,熟练使用相关API 理解各个API背后的相关算法原理,每个参数意义 有能力解决使用应用场景问题,大量工程代码经验分享 掌握图像处理与视频分析,图像分析与测量编码与开发技巧

    67人学习 贾志刚
    免费试看

http://www.ituring.com.cn/tupubarticle/2024

第 1 章 基本的图像操作和处理

本章讲解操作和处理图像的基础知识,将通过大量示例介绍处理图像所需的 Python 工具包,并介绍用于读取图像、图像转换和缩放、计算导数、画图和保存结果等的基本工具。这些工具的使用将贯穿本书的剩余章节。

1.1 PIL:Python图像处理类库

PIL(Python Imaging Library Python,图像处理类库)提供了通用的图像处理功能,以及大量有用的基本图像操作,比如图像缩放、裁剪、旋转、颜色转换等。PIL 是免费的,可以从 http://www.pythonware.com/products/pil/ 下载。

利用 PIL 中的函数,我们可以从大多数图像格式的文件中读取数据,然后写入最常见的图像格式文件中。PIL 中最重要的模块为 Image。要读取一幅图像,可以使用:

from PIL import Image

pil_im = Image.open('empire.jpg')

上述代码的返回值 pil_im 是一个 PIL 图像对象。

图像的颜色转换可以使用 convert() 方法来实现。要读取一幅图像,并将其转换成灰度图像,只需要加上 convert('L'),如下所示:

pil_im = Image.open('empire.jpg').convert('L')

在 PIL 文档中有一些例子,参见 http://www.pythonware.com/library/pil/handbook/index.htm。这些例子的输出结果如图 1-1 所示。

图 1-1:用 PIL 处理图像的例子

1.1.1 转换图像格式

通过 save() 方法,PIL 可以将图像保存成多种格式的文件。下面的例子从文件名列表(filelist)中读取所有的图像文件,并转换成 JPEG 格式:

from PIL import Image
import os

for infile in filelist:
  outfile = os.path.splitext(infile)[0] + ".jpg"
  if infile != outfile:
    try:
      Image.open(infile).save(outfile)
    except IOError:
      print "cannot convert", infile

PIL 的 open() 函数用于创建 PIL 图像对象,save() 方法用于保存图像到具有指定文件名的文件。除了后缀变为“.jpg”,上述代码的新文件名和原文件名相同。PIL 是个足够智能的类库,可以根据文件扩展名来判定图像的格式。PIL 函数会进行简单的检查,如果文件不是 JPEG 格式,会自动将其转换成 JPEG 格式;如果转换失败,它会在控制台输出一条报告失败的消息。

本书会处理大量图像列表。下面将创建一个包含文件夹中所有图像文件的文件名列表。首先新建一个文件,命名为 imtools.py,来存储一些经常使用的图像操作,然后将下面的函数添加进去:

import os
def get_imlist(path):

""" 返回目录中所有JPG 图像的文件名列表"""

return [os.path.join(path,f) for f in os.listdir(path) if f.endswith('.jpg')]

现在,回到 PIL。

1.1.2 创建缩略图

使用 PIL 可以很方便地创建图像的缩略图。thumbnail() 方法接受一个元组参数(该参数指定生成缩略图的大小),然后将图像转换成符合元组参数指定大小的缩略图。例如,创建最长边为 128 像素的缩略图,可以使用下列命令:

pil_im.thumbnail((128,128))

1.1.3 复制和粘贴图像区域

使用 crop() 方法可以从一幅图像中裁剪指定区域:

box = (100,100,400,400)
region = pil_im.crop(box)

该区域使用四元组来指定。四元组的坐标依次是(左,上,右,下)。PIL 中指定坐标系的左上角坐标为(0,0)。我们可以旋转上面代码中获取的区域,然后使用 paste() 方法将该区域放回去,具体实现如下:

region = region.transpose(Image.ROTATE_180)
pil_im.paste(region,box)

1.1.4 调整尺寸和旋转

要调整一幅图像的尺寸,我们可以调用 resize() 方法。该方法的参数是一个元组,用来指定新图像的大小:

out = pil_im.resize((128,128))

要旋转一幅图像,可以使用逆时针方式表示旋转角度,然后调用 rotate() 方法:

out = pil_im.rotate(45)

上述例子的输出结果如图 1-1 所示。最左端是原始图像,然后是灰度图像、粘贴有旋转后裁剪图像的原始图像,最后是缩略图。

1.2 Matplotlib

我们处理数学运算、绘制图表,或者在图像上绘制点、直线和曲线时,Matplotlib 是个很好的类库,具有比 PIL 更强大的绘图功能。Matplotlib 可以绘制出高质量的图表,就像本书中的许多插图一样。Matplotlib 中的 PyLab 接口包含很多方便用户创建图像的函数。Matplotlib 是开源工具,可以从 http://matplotlib.sourceforge.net/ 免费下载。该链接中包含非常详尽的使用说明和教程。下面的例子展示了本书中需要使用的大部分函数。

1.2.1 绘制图像、点和线

尽管 Matplotlib 可以绘制出较好的条形图、饼状图、散点图等,但是对于大多数计算机视觉应用来说,仅仅需要用到几个绘图命令。最重要的是,我们想用点和线来表示一些事物,比如兴趣点、对应点以及检测出的物体。下面是用几个点和一条线绘制图像的例子:

from PIL import Image
from pylab import *

# 读取图像到数组中
im = array(Image.open('empire.jpg'))

# 绘制图像
imshow(im)

# 一些点
x = [100,100,400,400]
y = [200,500,200,500]

# 使用红色星状标记绘制点
plot(x,y,'r*')

# 绘制连接前两个点的线
plot(x[:2],y[:2])

# 添加标题,显示绘制的图像
title('Plotting: "empire.jpg"')
show()

上面的代码首先绘制出原始图像,然后在 x 和 y 列表中给定点的 x 坐标和 y 坐标上绘制出红色星状标记点,最后在两个列表表示的前两个点之间绘制一条线段(默认为蓝色)。该例子的绘制结果如图 1-2 所示。show() 命令首先打开图形用户界面(GUI),然后新建一个图像窗口。该图形用户界面会循环阻断脚本,然后暂停,直到最后一个图像窗口关闭。在每个脚本里,你只能调用一次 show() 命令,而且通常是在脚本的结尾调用。注意,在 PyLab 库中,我们约定图像的左上角为坐标原点。

图像的坐标轴是一个很有用的调试工具;但是,如果你想绘制出较美观的图像,加上下列命令可以使坐标轴不显示:

axis('off')

上面的命令将绘制出如图 1-2 右边所示的图像。

图 1-2:Matplotlib 绘图示例。带有坐标轴和不带坐标轴的包含点和一条线段的图像

在绘图时,有很多选项可以控制图像的颜色和样式。最有用的一些短命令如表 1-1、表 1-2 和表 1-3 所示。使用方法见下面的例子:

plot(x,y)         # 默认为蓝色实线
plot(x,y,'r*')    # 红色星状标记
plot(x,y,'go-')   # 带有圆圈标记的绿线
plot(x,y,'ks:')   # 带有正方形标记的黑色虚线

表1-1:用PyLab库绘图的基本颜色格式命令

颜色

 

'b'

蓝色

'g'

绿色

'r'

红色

'c'

青色

'm'

品红

'y'

黄色

'k'

黑色

'w'

白色

表1-2:用PyLab库绘图的基本线型格式命令

线型

 

'-'

实线

'--'

虚线

':'

点线

表1-3:用PyLab库绘图的基本绘制标记格式命令

标记

 

'.'

'o'

圆圈

's'

正方形

'*'

星形

'+'

加号

'x'

叉号

1.2.2 图像轮廓和直方图

下面来看两个特别的绘图示例:图像的轮廓和直方图。绘制图像的轮廓(或者其他二维函数的等轮廓线)在工作中非常有用。因为绘制轮廓需要对每个坐标 [x, y] 的像素值施加同一个阈值,所以首先需要将图像灰度化:

from PIL import Image
from pylab import *

# 读取图像到数组中
im = array(Image.open('empire.jpg').convert('L'))

# 新建一个图像
figure()
# 不使用颜色信息
gray()
# 在原点的左上角显示轮廓图像
contour(im, origin='image')
axis('equal')
axis('off')

像之前的例子一样,这里用 PIL 的 convert() 方法将图像转换成灰度图像。

图像的直方图用来表征该图像像素值的分布情况。用一定数目的小区间(bin)来指定表征像素值的范围,每个小区间会得到落入该小区间表示范围的像素数目。该(灰度)图像的直方图可以使用 hist() 函数绘制:

figure()
hist(im.flatten(),128)
show()

hist() 函数的第二个参数指定小区间的数目。需要注意的是,因为 hist() 只接受一维数组作为输入,所以我们在绘制图像直方图之前,必须先对图像进行压平处理。flatten() 方法将任意数组按照行优先准则转换成一维数组。图 1-3 为等轮廓线和直方图图像。

图 1-3:用 Matplotlib 绘制图像等轮廓线和直方图

1.2.3 交互式标注

有时用户需要和某些应用交互,例如在一幅图像中标记一些点,或者标注一些训练数据。PyLab 库中的 ginput() 函数就可以实现交互式标注。下面是一个简短的例子:

from PIL import Image
from pylab import *

im = array(Image.open('empire.jpg'))
imshow(im)
print 'Please click 3 points'
x = ginput(3)
print 'you clicked:',x
show()

上面的脚本首先绘制一幅图像,然后等待用户在绘图窗口的图像区域点击三次。程序将这些点击的坐标 [x, y] 自动保存在 x 列表里。

1.3 NumPy

NumPyhttp://www.scipy.org/NumPy/)是非常有名的 Python 科学计算工具包,其中包含了大量有用的思想,比如数组对象(用来表示向量、矩阵、图像等)以及线性代数函数。NumPy 中的数组对象几乎贯穿用于本书的所有例子中 1 数组对象可以帮助你实现数组中重要的操作,比如矩阵乘积、转置、解方程系统、向量乘积和归一化,这为图像变形、对变化进行建模、图像分类、图像聚类等提供了基础。

1PyLab 实际上包含 NumPy 的一些内容,如数组类型。这也是我们能够在 1.2 节使用数组类型的原因。

NumPy 可以从 http://www.scipy.org/Download 免费下载,在线说明文档(http://docs.scipy.org/doc/numpy/)包含了你可能遇到的大多数问题的答案。关于 NumPy 的更多内容,请参考开源书籍 [24]。

1.3.1 图像数组表示

在先前的例子中,当载入图像时,我们通过调用 array() 方法将图像转换成 NumPy 的数组对象,但当时并没有进行详细介绍。NumPy 中的数组对象是多维的,可以用来表示向量、矩阵和图像。一个数组对象很像一个列表(或者是列表的列表),但是数组中所有的元素必须具有相同的数据类型。除非创建数组对象时指定数据类型,否则数据类型会按照数据的类型自动确定。

对于图像数据,下面的例子阐述了这一点:

im = array(Image.open('empire.jpg'))
print im.shape, im.dtype

im = array(Image.open('empire.jpg').convert('L'),'f')
print im.shape, im.dtype

控制台输出结果如下所示:

(800, 569, 3) uint8
(800, 569) float32

每行的第一个元组表示图像数组的大小(行、列、颜色通道),紧接着的字符串表示数组元素的数据类型。因为图像通常被编码成无符号八位整数(uint8),所以在第一种情况下,载入图像并将其转换到数组中,数组的数据类型为“uint8”。在第二种情况下,对图像进行灰度化处理,并且在创建数组时使用额外的参数“f”;该参数将数据类型转换为浮点型。关于更多数据类型选项,可以参考图书 [24]。注意,由于灰度图像没有颜色信息,所以在形状元组中,它只有两个数值。

数组中的元素可以使用下标访问。位于坐标 ij,以及颜色通道 k 的像素值可以像下面这样访问:

value = im[i,j,k]

多个数组元素可以使用数组切片方式访问。切片方式返回的是以指定间隔下标访问该数组的元素值。下面是有关灰度图像的一些例子:

im[i,:] = im[j,:]      # 将第 j 行的数值赋值给第 i 行
im[:,i] = 100          # 将第 i 列的所有数值设为100
im[:100,:50].sum()     # 计算前100 行、前 50 列所有数值的和
im[50:100,50:100]      # 50~100 行,50~100 列(不包括第 100 行和第 100 列)
im[i].mean()           # 第 i 行所有数值的平均值
im[:,-1]               # 最后一列
im[-2,:] (or im[-2])   # 倒数第二行

注意,示例仅仅使用一个下标访问数组。如果仅使用一个下标,则该下标为行下标。注意,在最后几个例子中,负数切片表示从最后一个元素逆向计数。我们将会频繁地使用切片技术访问像素值,这也是一个很重要的思想。

我们有很多操作和方法来处理数组对象。本书将在使用到的地方逐一介绍。你可以查阅在线文档或者开源图书 [24] 获取更多信息。

1.3.2 灰度变换

将图像读入 NumPy 数组对象后,我们可以对它们执行任意数学操作。一个简单的例子就是图像的灰度变换。考虑任意函数 f,它将 0...255 区间(或者 0...1 区间)映射到自身(意思是说,输出区间的范围和输入区间的范围相同)。下面是关于灰度变换的一些例子:

from PIL import Image
from numpy import *

im = array(Image.open('empire.jpg').convert('L'))

im2 = 255 - im # 对图像进行反相处理

im3 = (100.0/255) * im + 100 # 将图像像素值变换到100...200 区间

im4 = 255.0 * (im/255.0)**2 # 对图像像素值求平方后得到的图像

第一个例子将灰度图像进行反相处理;第二个例子将图像的像素值变换到 100...200 区间;第三个例子对图像使用二次函数变换,使较暗的像素值变得更小。图 1-4 为所使用的变换函数图像。图 1-5 是输出的图像结果。你可以使用下面的命令查看图像中的最小和最大像素值:

print int(im.min()), int(im.max())

图 1-4:灰度变换示例。三个例子中所使用函数的图像,其中虚线表示恒等变换

图 1-5:灰度变换。对图像应用图 1-4 中的函数:f(x)=255-x 对图像进行反相处理(左);f(x)=(100/255)x+100 对图像进行变换(中);f(x)=255(x/255)2 对图像做二次变换(右)

如果试着对上面例子查看最小值和最大值,可以得到下面的输出结果:

2 255
0 253
100 200
0 255

array() 变换的相反操作可以使用 PIL 的 fromarray() 函数完成:

pil_im = Image.fromarray(im)

如果你通过一些操作将“uint8”数据类型转换为其他数据类型,比如之前例子中的 im3 或者 im4,那么在创建 PIL 图像之前,需要将数据类型转换回来:

pil_im = Image.fromarray(uint8(im))

如果你并不十分确定输入数据的类型,安全起见,应该先转换回来。注意,NumPy 总是将数组数据类型转换成能够表示数据的“最低”数据类型。对浮点数做乘积或除法操作会使整数类型的数组变成浮点类型。

1.3.3 图像缩放

NumPy 的数组对象是我们处理图像和数据的主要工具。想要对图像进行缩放处理没有现成简单的方法。我们可以使用之前 PIL 对图像对象转换的操作,写一个简单的用于图像缩放的函数。把下面的函数添加到 imtool.py 文件里:

def imresize(im,sz):
  """ 使用PIL 对象重新定义图像数组的大小"""
  pil_im = Image.fromarray(uint8(im))

  return array(pil_im.resize(sz))

我们将会在接下来的内容中使用这个函数。

1.3.4 直方图均衡化

图像灰度变换中一个非常有用的例子就是直方图均衡化。直方图均衡化是指将一幅图像的灰度直方图变平,使变换后的图像中每个灰度值的分布概率都相同。在对图像做进一步处理之前,直方图均衡化通常是对图像灰度值进行归一化的一个非常好的方法,并且可以增强图像的对比度。

在这种情况下,直方图均衡化的变换函数是图像中像素值的累积分布函数(cumulative distribution function,简写为 cdf,将像素值的范围映射到目标范围的归一化操作)。

下面的函数是直方图均衡化的具体实现。将这个函数添加到 imtool.py 里:

def histeq(im,nbr_bins=256):
  """ 对一幅灰度图像进行直方图均衡化"""

  # 计算图像的直方图
  imhist,bins = histogram(im.flatten(),nbr_bins,normed=True)
  cdf = imhist.cumsum() # cumulative distribution function
  cdf = 255 * cdf / cdf[-1] # 归一化

  # 使用累积分布函数的线性插值,计算新的像素值
  im2 = interp(im.flatten(),bins[:-1],cdf)

  return im2.reshape(im.shape), cdf

该函数有两个输入参数,一个是灰度图像,一个是直方图中使用小区间的数目。函数返回直方图均衡化后的图像,以及用来做像素值映射的累积分布函数。注意,函数中使用到累积分布函数的最后一个元素(下标为 -1),目的是将其归一化到 0...1 范围。你可以像下面这样使用该函数:

from PIL import Image
from numpy import *

im = array(Image.open('AquaTermi_lowcontrast.jpg').convert('L'))
im2,cdf = imtools.histeq(im)

图 1-6 和图 1-7 为上面直方图均衡化例子的结果。上面一行显示的分别是直方图均衡化之前和之后的灰度直方图,以及累积概率分布函数映射图像。可以看到,直方图均衡化后图像的对比度增强了,原先图像灰色区域的细节变得清晰。

图 1-6:直方图均衡化示例。左侧为原始图像和直方图,中间图为灰度变换函数,右侧为直方图均衡化后的图像和相应直方图

图 1-7:直方图均衡化示例。左侧为原始图像和直方图,中间图为灰度变换函数,右侧为直方图均衡化后的图像和相应直方图

1.3.5 图像平均

图像平均操作是减少图像噪声的一种简单方式,通常用于艺术特效。我们可以简单地从图像列表中计算出一幅平均图像。假设所有的图像具有相同的大小,我们可以将这些图像简单地相加,然后除以图像的数目,来计算平均图像。下面的函数可以用于计算平均图像,将其添加到 imtool.py 文件里:

def compute_average(imlist):
  """ 计算图像列表的平均图像"""

  # 打开第一幅图像,将其存储在浮点型数组中
  averageim = array(Image.open(imlist[0]), 'f')

  for imname in imlist[1:]:
    try:
      averageim += array(Image.open(imname))
    except:
      print imname + '...skipped'
  averageim /= len(imlist)

  # 返回uint8 类型的平均图像
  return array(averageim, 'uint8')

该函数包括一些基本的异常处理技巧,可以自动跳过不能打开的图像。我们还可以使用 mean() 函数计算平均图像。mean() 函数需要将所有的图像堆积到一个数组中;也就是说,如果有很多图像,该处理方式需要占用很多内存。我们将会在下一节中使用该函数。

1.3.6 图像的主成分分析(PCA)

PCA(Principal Component Analysis,主成分分析)是一个非常有用的降维技巧。它可以在使用尽可能少维数的前提下,尽量多地保持训练数据的信息,在此意义上是一个最佳技巧。即使是一幅 100×100 像素的小灰度图像,也有 10 000 维,可以看成 10 000 维空间中的一个点。一兆像素的图像具有百万维。由于图像具有很高的维数,在许多计算机视觉应用中,我们经常使用降维操作。PCA 产生的投影矩阵可以被视为将原始坐标变换到现有的坐标系,坐标系中的各个坐标按照重要性递减排列。

为了对图像数据进行 PCA 变换,图像需要转换成一维向量表示。我们可以使用 NumPy 类库中的 flatten() 方法进行变换。

将变平的图像堆积起来,我们可以得到一个矩阵,矩阵的一行表示一幅图像。在计算主方向之前,所有的行图像按照平均图像进行了中心化。我们通常使用 SVD(Singular Value Decomposition,奇异值分解)方法来计算主成分;但当矩阵的维数很大时,SVD 的计算非常慢,所以此时通常不使用 SVD 分解。下面就是 PCA 操作的代码:

from PIL import Image
from numpy import *

def pca(X):
  """ 主成分分析:
    输入:矩阵X ,其中该矩阵中存储训练数据,每一行为一条训练数据
    返回:投影矩阵(按照维度的重要性排序)、方差和均值"""

  # 获取维数
  num_data,dim = X.shape

  # 数据中心化
  mean_X = X.mean(axis=0)
  X = X - mean_X

if dim>num_data:
  # PCA- 使用紧致技巧
  M = dot(X,X.T) # 协方差矩阵
  e,EV = linalg.eigh(M) # 特征值和特征向量
  tmp = dot(X.T,EV).T # 这就是紧致技巧
  V = tmp[::-1] # 由于最后的特征向量是我们所需要的,所以需要将其逆转
  S = sqrt(e)[::-1] # 由于特征值是按照递增顺序排列的,所以需要将其逆转
  for i in range(V.shape[1]):
    V[:,i] /= S
else:
  # PCA- 使用SVD 方法
  U,S,V = linalg.svd(X)
  V = V[:num_data] # 仅仅返回前nun_data 维的数据才合理

# 返回投影矩阵、方差和均值
return V,S,mean_X

该函数首先通过减去每一维的均值将数据中心化,然后计算协方差矩阵对应最大特征值的特征向量,此时可以使用简明的技巧或者 SVD 分解。这里我们使用了 range() 函数,该函数的输入参数为一个整数 n,函数返回整数 0...(n-1) 的一个列表。你也可以使用 arange() 函数来返回一个数组,或者使用 xrange() 函数返回一个产生器(可能会提升速度)。我们在本书中贯穿使用 range() 函数。

如果数据个数小于向量的维数,我们不用 SVD 分解,而是计算维数更小的协方差矩阵 XXT 的特征向量。通过仅计算对应前 kk 是降维后的维数)最大特征值的特征向量,可以使上面的 PCA 操作更快。由于篇幅所限,有兴趣的读者可以自行探索。矩阵 V 的每行向量都是正交的,并且包含了训练数据方差依次减少的坐标方向。

我们接下来对字体图像进行 PCA 变换。fontimages.zip 文件包含采用不同字体的字符 a 的缩略图。所有的 2359 种字体可以免费下载 2。假定这些图像的名称保存在列表 imlist 中,跟之前的代码一起保存传在 pca.py 文件中,我们可以使用下面的脚本计算图像的主成分:

2免费字体图像库由 Martin Solli 收集并上传(http://webstaff.itn.liu.se/~marso/)。

from PIL import Image
from numpy import *
from pylab import *
import pca

im = array(Image.open(imlist[0])) # 打开一幅图像,获取其大小
m,n = im.shape[0:2] # 获取图像的大小
imnbr = len(imlist) # 获取图像的数目

# 创建矩阵,保存所有压平后的图像数据
immatrix = array([array(Image.open(im)).flatten()
               for im in imlist],'f')

# 执行 PCA 操作
V,S,immean = pca.pca(immatrix)

# 显示一些图像(均值图像和前 7 个模式)
figure()
gray()
subplot(2,4,1)
imshow(immean.reshape(m,n))
for i in range(7):
  subplot(2,4,i+2)
  imshow(V[i].reshape(m,n))

show()

注意,图像需要从一维表示重新转换成二维图像;可以使用 reshape() 函数。如图 1-8 所示,运行该例子会在一个绘图窗口中显示 8 个图像。这里我们使用了 PyLab 库的 subplot() 函数在一个窗口中放置多个图像。

图 1-8:平均图像(左上)和前 7 个模式(具有最大方差的方向模式)

1.3.7 使用pickle模块

如果想要保存一些结果或者数据以方便后续使用,Python 中的 pickle 模块非常有用。pickle 模块可以接受几乎所有的 Python 对象,并且将其转换成字符串表示,该过程叫做封装(pickling)。从字符串表示中重构该对象,称为拆封(unpickling)。这些字符串表示可以方便地存储和传输。

我们来看一个例子。假设想要保存上一节字体图像的平均图像和主成分,可以这样来完成:

# 保存均值和主成分数据
f = open('font_pca_modes.pkl', 'wb')
pickle.dump(immean,f)
pickle.dump(V,f)
f.close()

在上述例子中,许多对象可以保存到同一个文件中。pickle 模块中有很多不同的协议可以生成 .pkl 文件;如果不确定的话,最好以二进制文件的形式读取和写入。在其他 Python 会话中载入数据,只需要如下使用 load() 方法:

# 载入均值和主成分数据
f = open('font_pca_modes.pkl', 'rb')
immean = pickle.load(f)
V = pickle.load(f)
f.close()

注意,载入对象的顺序必须和先前保存的一样。Python 中有个用 C 语言写的优化版本,叫做 cpickle 模块,该模块和标准 pickle 模块完全兼容。关于 pickle 模块的更多内容,参见 pickle 模块文档页 http://docs.python.org/library/pickle.html

在本书接下来的章节中,我们将使用 with 语句处理文件的读写操作。这是 Python 2.5 引入的思想,可以自动打开和关闭文件(即使在文件打开时发生错误)。下面的例子使用 with() 来实现保存和载入操作:

# 打开文件并保存
with open('font_pca_modes.pkl', 'wb') as f:
  pickle.dump(immean,f)
  pickle.dump(V,f)

# 打开文件并载入
with open('font_pca_modes.pkl', 'rb') as f:
  immean = pickle.load(f)
  V = pickle.load(f)

上面的例子乍看起来可能很奇怪,但 with() 确实是个很有用的思想。如果你不喜欢它,可以使用之前的 openclose 函数。

作为 pickle 的一种替代方式,NumPy 具有读写文本文件的简单函数。如果数据中不包含复杂的数据结构,比如在一幅图像上点击的点列表,NumPy 的读写函数会很有用。保存一个数组 x 到文件中,可以使用:

savetxt('test.txt',x,'%i')

最后一个参数表示应该使用整数格式。类似地,读取可以使用:

x = loadtxt('test.txt')

你可以从在线文档 http://docs.scipy.org/doc/numpy/reference/generated/numpy.loadtxt.html 了解更多内容。

最后,NumPy 有专门用于保存和载入数组的函数。你可以在上面的在线文档里查看关于 save() load() 的更多内容。

1.4 SciPy

SciPyhttp://scipy.org/) 是建立在 NumPy 基础上,用于数值运算的开源工具包。SciPy 提供很多高效的操作,可以实现数值积分、优化、统计、信号处理,以及对我们来说最重要的图像处理功能。接下来,本节会介绍 SciPy 中大量有用的模块。SciPy 是个开源工具包,可以从 http://scipy.org/Download 下载。

1.4.1 图像模糊

图像的高斯模糊是非常经典的图像卷积例子。本质上,图像模糊就是将(灰度)图像 I 和一个高斯核进行卷积操作:

Iσ = I*Gσ

其中 * 表示卷积操作;Gσ 是标准差为 σ 的二维高斯核,定义为 :

G_\sigma=\frac{1}{2\pi\sigma}e^{-(x^2+y^2)/2\sigma^2}

高斯模糊通常是其他图像处理操作的一部分,比如图像插值操作、兴趣点计算以及很多其他应用。

SciPy 有用来做滤波操作的 scipy.ndimage.filters 模块。该模块使用快速一维分离的方式来计算卷积。你可以像下面这样来使用它:

from PIL import Image
from numpy import *
from scipy.ndimage import filters

im = array(Image.open('empire.jpg').convert('L'))
im2 = filters.gaussian_filter(im,5)

上面 guassian_filter() 函数的最后一个参数表示标准差。

图 1-9 显示了随着 σ 的增加,一幅图像被模糊的程度。σ 越大,处理后的图像细节丢失越多。如果打算模糊一幅彩色图像,只需简单地对每一个颜色通道进行高斯模糊:

im = array(Image.open('empire.jpg'))
im2 = zeros(im.shape)
for i in range(3):
  im2[:,:,i] = filters.gaussian_filter(im[:,:,i],5)
im2 = uint8(im2)

在上面的脚本中,最后并不总是需要将图像转换成 uint8 格式,这里只是将像素值用八位来表示。我们也可以使用:

im2 = array(im2,'uint8')

来完成转换。

关于该模块更多的内容以及不同参数的选择,请查看 http://docs.scipy.org/doc/scipy/reference/ndimage.htmlSciPy 文档中的 scipy.ndimage 部分。

图 1-9:使用 scipy.ndimage.filters 模块进行高斯模糊:(a)原始灰度图像;(b)使用 σ=2 的高斯滤波器;(c)使用 σ=5 的高斯滤波器;(d)使用 σ=10 的高斯滤波器

1.4.2 图像导数

整本书中可以看到,在很多应用中图像强度的变化情况是非常重要的信息。强度的变化可以用灰度图像 I(对于彩色图像,通常对每个颜色通道分别计算导数)的 xy 方向导数 Ix Iy 进行描述。

图像的梯度向量为∇I = [Ix, Iy]T。梯度有两个重要的属性,一是梯度的大小

\left|\boldsymbol{\nabla I}\right|=\sqrt{{\boldsymbol{I}_x}^2+{\boldsymbol{I}_y}^2}

它描述了图像强度变化的强弱,一是梯度的角度

α=arctan2(Iy, Ix)

描述了图像中在每个点(像素)上强度变化最大的方向。NumPy 中的 arctan2() 函数返回弧度表示的有符号角度,角度的变化区间为 -π...π。

我们可以用离散近似的方式来计算图像的导数。图像导数大多数可以通过卷积简单地实现:

Ix=I*DxIy=I*Dy

对于 DxDy,通常选择 Prewitt 滤波器:

D_x=\begin{vmatrix}-1&0&1\\-1&0&1\\ -1&0&1\end{vmatrix}D_y=\begin{vmatrix}-1&-1&-1\\0&0&0\\ 1&1&1\end{vmatrix}

或者 Sobel 滤波器:

D_x=\begin{vmatrix}-1&0&1\\-2&0&2\\ -1&0&1\end{vmatrix}D_y=\begin{vmatrix}-1&-2&-1\\0&0&0\\ 1&2&1\end{vmatrix}

这些导数滤波器可以使用 scipy.ndimage.filters 模块的标准卷积操作来简单地实现,例如:

from PIL import Image
from numpy import *
from scipy.ndimage import filters

im = array(Image.open('empire.jpg').convert('L'))

# Sobel 导数滤波器
imx = zeros(im.shape)
filters.sobel(im,1,imx)

imy = zeros(im.shape)
filters.sobel(im,0,imy)

magnitude = sqrt(imx**2+imy**2)

上面的脚本使用 Sobel 滤波器来计算 xy 的方向导数,以及梯度大小。sobel() 函数的第二个参数表示选择 x 或者 y 方向导数,第三个参数保存输出的变量。图 1-10 显示了用 Sobel 滤波器计算出的导数图像。在两个导数图像中,正导数显示为亮的像素,负导数显示为暗的像素。灰色区域表示导数的值接近于零。

图 1-10:使用 Sobel 导数滤波器计算导数图像:(a)原始灰度图像;(b)x 导数图像;(c)y 导数图像;(d)梯度大小图像

上述计算图像导数的方法有一些缺陷:在该方法中,滤波器的尺度需要随着图像分辨率的变化而变化。为了在图像噪声方面更稳健,以及在任意尺度上计算导数,我们可以使用高斯导数滤波器:

Ix=I*GσxIy=I*Gσy

其中,GσxGσy 表示 Gσxy 方向上的导数,Gσ 为标准差为 σ 的高斯函数。

我们之前用于模糊的 filters.gaussian_filter() 函数可以接受额外的参数,用来计算高斯导数。可以简单地按照下面的方式来处理:

sigma = 5 # 标准差

imx = zeros(im.shape)
filters.gaussian_filter(im, (sigma,sigma), (0,1), imx)

imy = zeros(im.shape)
filters.gaussian_filter(im, (sigma,sigma), (1,0), imy)

该函数的第三个参数指定对每个方向计算哪种类型的导数,第二个参数为使用的标准差。你可以查看相应文档了解详情。图 1-11 显示了不同尺度下的导数图像和梯度大小。你可以和图 1-9 中做相同尺度模糊的图像做比较。

图 1-11:使用高斯导数计算图像导数:x 导数图像(上),y 导数图像(中),以及梯度大小图像(下);(a)为原始灰度图像,(b)为使用 σ=2 的高斯导数滤波器处理后的图像,(c)为使 用 σ=5 的高斯导数滤波器处理后的图像,(d)为使用 σ=10 的高斯导数滤波器处理后的图像

1.4.3 形态学:对象计数

形态学(或数学形态学)是度量和分析基本形状的图像处理方法的基本框架与集合。形态学通常用于处理二值图像,但是也能够用于灰度图像。二值图像是指图像的每个像素只能取两个值,通常是 0 和 1。二值图像通常是,在计算物体的数目,或者度量其大小时,对一幅图像进行阈值化后的结果。你可以从 http://en.wikipedia.org/wiki/Mathematical_morphology 大体了解形态学及其处理图像的方式。

scipy.ndimage 中的 morphology 模块可以实现形态学操作。你可以使用 scipy.ndimage 中的 measurements 模块来实现二值图像的计数和度量功能。下面通过一个简单的例子介绍如何使用它们。

考虑在图 1-12a3 里的二值图像,计算该图像中的对象个数可以通过下面的脚本实现:

3这个图像实际上是图像“分割”后的结果。如果你想知道该图像是如何创建的,可以查看 9.3 节。

from scipy.ndimage import measurements,morphology

# 载入图像,然后使用阈值化操作,以保证处理的图像为二值图像
im = array(Image.open('houses.png').convert('L'))
im = 1*(im<128)

labels, nbr_objects = measurements.label(im)
print "Number of objects:", nbr_objects

上面的脚本首先载入该图像,通过阈值化方式来确保该图像是二值图像。通过和 1 相乘,脚本将布尔数组转换成二进制表示。然后,我们使用 label() 函数寻找单个的物体,并且按照它们属于哪个对象将整数标签给像素赋值。图 1-12b 是 labels 数组的图像。图像的灰度值表示对象的标签。可以看到,在一些对象之间有一些小的连接。进行二进制开(binary open)操作,我们可以将其移除:

# 形态学开操作更好地分离各个对象
im_open = morphology.binary_opening(im,ones((9,5)),iterations=2)

labels_open, nbr_objects_open = measurements.label(im_open)
print "Number of objects:", nbr_objects_open

binary_opening() 函数的第二个参数指定一个数组结构元素。该数组表示以一个像素为中心时,使用哪些相邻像素。在这种情况下,我们在 y 方向上使用 9 个像素(上面 4 个像素、像素本身、下面 4 个像素),在 x 方向上使用 5 个像素。你可以指定任意数组为结构元素,数组中的非零元素决定使用哪些相邻像素。参数 iterations 决定执行该操作的次数。你可以尝试使用不同的迭代次数 iterations 值,看一下对象的数目如何变化。你可以在图 1-12c 与图 1-12d 中查看经过开操作后的图像,以及相应的标签图像。正如你想象的一样,binary_closing() 函数实现相反的操作。我们将该函数和在 morphologymeasurements 模块中的其他函数的用法留作练习。你可以从 scipy.ndimage 模块文档 http://docs.scipy.org/doc/scipy/reference/ndimage.html 中了解关于这些函数的更多知识。

图 1-12:形态学示例。使用二值开操作将对象分开,然后计算物体的数目:(a)为原始二值图像;(b)为对应原始图像的标签图像,其中灰度值表示物体的标签;(c)为使用开操作后的二值图像;(d)为开操作后图像的标签图像

1.4.4 一些有用的SciPy模块

SciPy 中包含一些用于输入和输出的实用模块。下面介绍其中两个模块:iomisc

  1. 读写.mat文件

    如果你有一些数据,或者在网上下载到一些有趣的数据集,这些数据以 Matlab 的 .mat 文件格式存储,那么可以使用 scipy.io 模块进行读取。

    data = scipy.io.loadmat('test.mat')
    

    上面代码中,data 对象包含一个字典,字典中的键对应于保存在原始 .mat 文件中的变量名。由于这些变量是数组格式的,因此可以很方便地保存到 .mat 文件中。你仅需创建一个字典(其中要包含你想要保存的所有变量),然后使用 savemat() 函数:

    data = {}
    data['x'] = x
    scipy.io.savemat('test.mat',data)
    

    因为上面的脚本保存的是数组 x,所以当读入到 Matlab 中时,变量的名字仍为 x。关于 scipy.io 模块的更多内容,请参见在线文档 http://docs.scipy.org/doc/scipy/reference/io.html

  2. 以图像形式保存数组

    因为我们需要对图像进行操作,并且需要使用数组对象来做运算,所以将数组直接保存为图像文件 4 非常有用。本书中的很多图像都是这样的创建的。

    imsave() 函数可以从 scipy.misc 模块中载入。要将数组 im 保存到文件中,可以使用下面的命令:

    from scipy.misc import imsave
    imsave('test.jpg',im)
    

    scipy.misc 模块同样包含了著名的 Lena 测试图像:

    lena = scipy.misc.lena()
    

    该脚本返回一个 512×512 的灰度图像数组。

4所有 Pylab 图均可保存为多种图像格式,方法是点击图像窗口中的“保存”按钮。

1.5 高级示例:图像去噪

我们通过一个非常实用的例子——图像的去噪——来结束本章。图像去噪是在去除图像噪声的同时,尽可能地保留图像细节和结构的处理技术。我们这里使用 ROF(Rudin-Osher-Fatemi)去噪模型。该模型最早出现在文献 [28] 中。图像去噪对于很多应用来说都非常重要;这些应用范围很广,小到让你的假期照片看起来更漂亮,大到提高卫星图像的质量。ROF 模型具有很好的性质:使处理后的图像更平滑,同时保持图像边缘和结构信息。

ROF 模型的数学基础和处理技巧非常高深,不在本书讲述范围之内。在讲述如何基于 Chambolle 提出的算法 [5] 实现 ROF 求解器之前,本书首先简要介绍一下 ROF 模型。

一幅(灰度)图像 I全变差(Total Variation,TV)定义为梯度范数之和。在连续表示的情况下,全变差表示为:

J(\boldsymbol{I})=\int\left|\nabla\boldsymbol{I}\right|\text{dx}            (1.1)

在离散表示的情况下,全变差表示为:

J(\boldsymbol{I})=\sum_{\text{x}}\left|\nabla\boldsymbol{I}\right|

其中,上面的式子是在所有图像坐标 x=[x, y] 上取和。

在 Chambolle 提出的 ROF 模型里,目标函数为寻找降噪后的图像 U,使下式最小:

\min_U\left|\left|\boldsymbol{I}-\boldsymbol{U}\right|\right|^2+2\lambda J(\boldsymbol{U}),

其中范数 ||I-U|| 是去噪后图像 U 和原始图像 I 差异的度量。也就是说,本质上该模型使去噪后的图像像素值“平坦”变化,但是在图像区域的边缘上,允许去噪后的图像像素值“跳跃”变化。

按照论文 [5] 中的算法,我们可以按照下面的代码实现 ROF 模型去噪:

from numpy import *

def denoise(im,U_init,tolerance=0.1,tau=0.125,tv_weight=100):
  """ 使用A. Chambolle(2005)在公式(11)中的计算步骤实现Rudin-Osher-Fatemi(ROF)去噪模型

    输入:含有噪声的输入图像(灰度图像)、U 的初始值、TV 正则项权值、步长、停业条件

    输出:去噪和去除纹理后的图像、纹理残留"""

m,n = im.shape # 噪声图像的大小

# 初始化
U = U_init
Px = im # 对偶域的x 分量
Py = im # 对偶域的y 分量
error = 1

while (error > tolerance):
  Uold = U

  # 原始变量的梯度
  GradUx = roll(U,-1,axis=1)-U # 变量U 梯度的x 分量
  GradUy = roll(U,-1,axis=0)-U # 变量U 梯度的y 分量

  # 更新对偶变量
  PxNew = Px + (tau/tv_weight)*GradUx
  PyNew = Py + (tau/tv_weight)*GradUy
  NormNew = maximum(1,sqrt(PxNew**2+PyNew**2))

  Px = PxNew/NormNew # 更新x 分量(对偶)
  Py = PyNew/NormNew # 更新y 分量(对偶)

  # 更新原始变量
  RxPx = roll(Px,1,axis=1) # 对x 分量进行向右x 轴平移
  RyPy = roll(Py,1,axis=0) # 对y 分量进行向右y 轴平移

  DivP = (Px-RxPx)+(Py-RyPy) # 对偶域的散度
  U = im + tv_weight*DivP # 更新原始变量

  # 更新误差
  error = linalg.norm(U-Uold)/sqrt(n*m);

return U,im-U # 去噪后的图像和纹理残余

在这个例子中,我们使用了 roll() 函数。顾名思义,在一个坐标轴上,它循环“滚动”数组中的元素值。该函数可以非常方便地计算邻域元素的差异,比如这里的导数。我们还使用了 linalg.norm() 函数,该函数可以衡量两个数组间(这个例子中是指图像矩阵 UUold)的差异。我们将这个 denoise() 函数保存到 rof.py 文件中。

下面使用一个合成的噪声图像示例来说明如何使用该函数:

from numpy import *
from numpy import random
from scipy.ndimage import filters
import rof

# 使用噪声创建合成图像
im = zeros((500,500))
im[100:400,100:400] = 128
im[200:300,200:300] = 255
im = im + 30*random.standard_normal((500,500))

U,T = rof.denoise(im,im)
G = filters.gaussian_filter(im,10)

# 保存生成结果
from scipy.misc import imsave
imsave('synth_rof.pdf',U)
imsave('synth_gaussian.pdf',G)

原始图像和图像的去噪结果如图 1-13 所示。正如你所看到的,ROF 算法去噪后的图像很好地保留了图像的边缘信息。

图 1-13:使用 ROF 模型对合成图像去噪:(a)为原始噪声图像;(b)为经过高斯模糊的图像(σ=10);(c)为经过 ROF 模型去噪后的图像

下面看一下在实际图像中使用 ROF 模型去噪的效果:

from PIL import Image
from pylab import *
import rof

im = array(Image.open('empire.jpg').convert('L'))
U,T = rof.denoise(im,im)

figure()
gray()
imshow(U)
axis('equal')
axis('off')
show()

经过 ROF 去噪后的图像如图 1-14c 所示。为了方便比较,该图中同样显示了模糊后的图像。可以看到,ROF 去噪后的图像保留了边缘和图像的结构信息,同时模糊了“噪声”。

图 1-14:使用 ROF 模型对灰度图像去噪:(a)为原始噪声图像;(b)为经过高斯模糊的图像(σ=5);(c)为经过 ROF 模型去噪后的图像

练习

  1. 如图 1-9 所示,将一幅图像进行高斯模糊处理。随着 σ 的增加,绘制出图像轮廓。在你绘制出的图中,图像的轮廓有何变化?你能解释为什么会发生这些变化吗?

  2. 通过将图像模糊化,然后从原始图像中减去模糊图像,来实现反锐化图像掩模操作(http://en.wikipedia.org/wiki/Unsharp_masking)。反锐化图像掩模操作可以实现图像锐化效果。试在彩色和灰度图像上使用反锐化图像掩模操作,观察该操作的效果。

  3. 除了直方图均衡化,商图像是另一种图像归一化的方法。商图像可以通过除以模糊后的图像 I/(I* Gσ) 获得。尝试使用该方法,并使用一些样本图像进行验证。

  4. 使用图像梯度,编写一个在图像中获得简单物体(例如,白色背景中的正方形)轮廓的函数。

  5. 使用梯度方向和大小检测图像中的线段。估计线段的长度以及线段的参数,并在原始图像中重新绘制该线段。

  6. 使用 label() 函数处理二值化图像,并使用直方图和标签图像绘制图像中物体的大小分布。

  7. 使用形态学操作处理阈值化图像。在发现一些参数能够产生好的结果后,使用 morphology 模块里面的 center_of_mass() 函数寻找每个物体的中心坐标,将其在图像中绘制出来。

代码示例约定

从第 2 章起,我们假定 PIL、NumPyMatplotlib 都包括在你所创建的每个文件和每个代码例子的开头:

from PIL import Image
from numpy import *
from pylab import *

这种约定使得示例代码更清晰,同时也便于读者理解。除此之外,我们使用 SciPy 模块时,将会在代码示例中显式声明。

一些纯化论者会反对这种将全体模块导入的方式,坚持如下使用方式:

import numpy as np
import matplotlib.pyplot as plt

这种方式能够保持命名空间(知道每个函数从哪儿来)。因为我们不需要 PyLab 中的 NumPy 部分,所以该例子只从 Matplotlib 中导入 pyplot 部分。纯化论者和经验丰富的程序员们知道这些区别,他们能够选择自己喜欢的方式。但是,为了使本书的内容和例子更容易被读者接受,我们不打算这样做。

请读者注意。


2017-02-10 15:28:58 xiao_lxl 阅读数 10040
  • 直线类型与线段

    掌握OpenCV核心模块,熟练使用相关API 理解各个API背后的相关算法原理,每个参数意义 有能力解决使用应用场景问题,大量工程代码经验分享 掌握图像处理与视频分析,图像分析与测量编码与开发技巧

    67人学习 贾志刚
    免费试看


转自 http://www.ituring.com.cn/tupubarticle/2024


第 1 章 基本的图像操作和处理

本章讲解操作和处理图像的基础知识,将通过大量示例介绍处理图像所需的 Python 工具包,并介绍用于读取图像、图像转换和缩放、计算导数、画图和保存结果等的基本工具。这些工具的使用将贯穿本书的剩余章节。

1.1 PIL:Python图像处理类库

PIL(Python Imaging Library,图像处理类库)提供了通用的图像处理功能,以及大量有用的基本图像操作,比如图像缩放、裁剪、旋转、颜色转换等。PIL 是免费的,可以从http://www.pythonware.com/products/pil/ 下载。

利用 PIL 中的函数,我们可以从大多数图像格式的文件中读取数据,然后写入最常见的图像格式文件中。PIL 中最重要的模块为 Image。要读取一幅图像,可以使用:

from PIL import Image

pil_im = Image.open('empire.jpg')

上述代码的返回值 pil_im 是一个 PIL 图像对象。

图像的颜色转换可以使用 convert() 方法来实现。要读取一幅图像,并将其转换成灰度图像,只需要加上 convert('L'),如下所示:

pil_im = Image.open('empire.jpg').convert('L')

在 PIL 文档中有一些例子,参见http://www.pythonware.com/library/pil/handbook/index.htm。这些例子的输出结果如图 1-1 所示。

图 1-1:用 PIL 处理图像的例子

1.1.1 转换图像格式

通过 save() 方法,PIL 可以将图像保存成多种格式的文件。下面的例子从文件名列表(filelist)中读取所有的图像文件,并转换成 JPEG 格式:

from PIL import Image
import os

for infile in filelist:
  outfile = os.path.splitext(infile)[0] + ".jpg"
  if infile != outfile:
    try:
      Image.open(infile).save(outfile)
    except IOError:
      print "cannot convert", infile

PIL 的 open() 函数用于创建 PIL 图像对象,save() 方法用于保存图像到具有指定文件名的文件。除了后缀变为“.jpg”,上述代码的新文件名和原文件名相同。PIL 是个足够智能的类库,可以根据文件扩展名来判定图像的格式。PIL 函数会进行简单的检查,如果文件不是 JPEG 格式,会自动将其转换成 JPEG 格式;如果转换失败,它会在控制台输出一条报告失败的消息。

本书会处理大量图像列表。下面将创建一个包含文件夹中所有图像文件的文件名列表。首先新建一个文件,命名为 imtools.py,来存储一些经常使用的图像操作,然后将下面的函数添加进去:

import os
def get_imlist(path):

""" 返回目录中所有JPG 图像的文件名列表"""

return [os.path.join(path,f) for f in os.listdir(path) if f.endswith('.jpg')]

现在,回到 PIL。

1.1.2 创建缩略图

使用 PIL 可以很方便地创建图像的缩略图。thumbnail() 方法接受一个元组参数(该参数指定生成缩略图的大小),然后将图像转换成符合元组参数指定大小的缩略图。例如,创建最长边为 128 像素的缩略图,可以使用下列命令:

pil_im.thumbnail((128,128))

1.1.3 复制和粘贴图像区域

使用 crop() 方法可以从一幅图像中裁剪指定区域:

box = (100,100,400,400)
region = pil_im.crop(box)

该区域使用四元组来指定。四元组的坐标依次是(左,上,右,下)。PIL 中指定坐标系的左上角坐标为(0,0)。我们可以旋转上面代码中获取的区域,然后使用 paste() 方法将该区域放回去,具体实现如下:

region = region.transpose(Image.ROTATE_180)
pil_im.paste(region,box)

1.1.4 调整尺寸和旋转

要调整一幅图像的尺寸,我们可以调用 resize() 方法。该方法的参数是一个元组,用来指定新图像的大小:

out = pil_im.resize((128,128))

要旋转一幅图像,可以使用逆时针方式表示旋转角度,然后调用 rotate() 方法:

out = pil_im.rotate(45)

上述例子的输出结果如图 1-1 所示。最左端是原始图像,然后是灰度图像、粘贴有旋转后裁剪图像的原始图像,最后是缩略图。

1.2 Matplotlib

我们处理数学运算、绘制图表,或者在图像上绘制点、直线和曲线时,Matplotlib 是个很好的类库,具有比 PIL 更强大的绘图功能。Matplotlib 可以绘制出高质量的图表,就像本书中的许多插图一样。Matplotlib 中的 PyLab 接口包含很多方便用户创建图像的函数。Matplotlib 是开源工具,可以从 http://matplotlib.sourceforge.net/ 免费下载。该链接中包含非常详尽的使用说明和教程。下面的例子展示了本书中需要使用的大部分函数。

1.2.1 绘制图像、点和线

尽管 Matplotlib 可以绘制出较好的条形图、饼状图、散点图等,但是对于大多数计算机视觉应用来说,仅仅需要用到几个绘图命令。最重要的是,我们想用点和线来表示一些事物,比如兴趣点、对应点以及检测出的物体。下面是用几个点和一条线绘制图像的例子:

from PIL import Image
from pylab import *

# 读取图像到数组中
im = array(Image.open('empire.jpg'))

# 绘制图像
imshow(im)

# 一些点
x = [100,100,400,400]
y = [200,500,200,500]

# 使用红色星状标记绘制点
plot(x,y,'r*')

# 绘制连接前两个点的线
plot(x[:2],y[:2])

# 添加标题,显示绘制的图像
title('Plotting: "empire.jpg"')
show()

上面的代码首先绘制出原始图像,然后在 x 和 y 列表中给定点的 x 坐标和 y 坐标上绘制出红色星状标记点,最后在两个列表表示的前两个点之间绘制一条线段(默认为蓝色)。该例子的绘制结果如图 1-2 所示。show() 命令首先打开图形用户界面(GUI),然后新建一个图像窗口。该图形用户界面会循环阻断脚本,然后暂停,直到最后一个图像窗口关闭。在每个脚本里,你只能调用一次show() 命令,而且通常是在脚本的结尾调用。注意,在 PyLab 库中,我们约定图像的左上角为坐标原点。

图像的坐标轴是一个很有用的调试工具;但是,如果你想绘制出较美观的图像,加上下列命令可以使坐标轴不显示:

axis('off')

上面的命令将绘制出如图 1-2 右边所示的图像。

图 1-2:Matplotlib 绘图示例。带有坐标轴和不带坐标轴的包含点和一条线段的图像

在绘图时,有很多选项可以控制图像的颜色和样式。最有用的一些短命令如表 1-1、表 1-2 和表 1-3 所示。使用方法见下面的例子:

plot(x,y)         # 默认为蓝色实线
plot(x,y,'r*')    # 红色星状标记
plot(x,y,'go-')   # 带有圆圈标记的绿线
plot(x,y,'ks:')   # 带有正方形标记的黑色点线

表1-1:用PyLab库绘图的基本颜色格式命令

颜色

 

'b'

蓝色

'g'

绿色

'r'

红色

'c'

青色

'm'

品红

'y'

黄色

'k'

黑色

'w'

白色

表1-2:用PyLab库绘图的基本线型格式命令

线型

 

'-'

实线

'--'

虚线

':'

点线

表1-3:用PyLab库绘图的基本绘制标记格式命令

标记

 

'.'

'o'

圆圈

's'

正方形

'*'

星形

'+'

加号

'x'

叉号

1.2.2 图像轮廓和直方图

下面来看两个特别的绘图示例:图像的轮廓和直方图。绘制图像的轮廓(或者其他二维函数的等轮廓线)在工作中非常有用。因为绘制轮廓需要对每个坐标 [x, y] 的像素值施加同一个阈值,所以首先需要将图像灰度化:

from PIL import Image
from pylab import *

# 读取图像到数组中
im = array(Image.open('empire.jpg').convert('L'))

# 新建一个图像
figure()
# 不使用颜色信息
gray()
# 在原点的左上角显示轮廓图像
contour(im, origin='image')
axis('equal')
axis('off')

像之前的例子一样,这里用 PIL 的 convert() 方法将图像转换成灰度图像。

图像的直方图用来表征该图像像素值的分布情况。用一定数目的小区间(bin)来指定表征像素值的范围,每个小区间会得到落入该小区间表示范围的像素数目。该(灰度)图像的直方图可以使用hist() 函数绘制:

figure()
hist(im.flatten(),128)
show()

hist() 函数的第二个参数指定小区间的数目。需要注意的是,因为 hist() 只接受一维数组作为输入,所以我们在绘制图像直方图之前,必须先对图像进行压平处理。flatten() 方法将任意数组按照行优先准则转换成一维数组。图 1-3 为等轮廓线和直方图图像。

图 1-3:用 Matplotlib 绘制图像等轮廓线和直方图

1.2.3 交互式标注

有时用户需要和某些应用交互,例如在一幅图像中标记一些点,或者标注一些训练数据。PyLab库中的 ginput() 函数就可以实现交互式标注。下面是一个简短的例子:

from PIL import Image
from pylab import *

im = array(Image.open('empire.jpg'))
imshow(im)
print 'Please click 3 points'
x = ginput(3)
print 'you clicked:',x
show()

上面的脚本首先绘制一幅图像,然后等待用户在绘图窗口的图像区域点击三次。程序将这些点击的坐标 [x, y] 自动保存在 x 列表里。

1.3 NumPy

NumPyhttp://www.scipy.org/NumPy/)是非常有名的 Python 科学计算工具包,其中包含了大量有用的思想,比如数组对象(用来表示向量、矩阵、图像等)以及线性代数函数。NumPy 中的数组对象几乎贯穿用于本书的所有例子中 1 数组对象可以帮助你实现数组中重要的操作,比如矩阵乘积、转置、解方程系统、向量乘积和归一化,这为图像变形、对变化进行建模、图像分类、图像聚类等提供了基础。

1PyLab 实际上包含 NumPy 的一些内容,如数组类型。这也是我们能够在 1.2 节使用数组类型的原因。

NumPy 可以从 http://www.scipy.org/Download 免费下载,在线说明文档(http://docs.scipy.org/doc/numpy/)包含了你可能遇到的大多数问题的答案。关于 NumPy 的更多内容,请参考开源书籍 [24]。

1.3.1 图像数组表示

在先前的例子中,当载入图像时,我们通过调用 array() 方法将图像转换成 NumPy 的数组对象,但当时并没有进行详细介绍。NumPy 中的数组对象是多维的,可以用来表示向量、矩阵和图像。一个数组对象很像一个列表(或者是列表的列表),但是数组中所有的元素必须具有相同的数据类型。除非创建数组对象时指定数据类型,否则数据类型会按照数据的类型自动确定。

对于图像数据,下面的例子阐述了这一点:

im = array(Image.open('empire.jpg'))
print im.shape, im.dtype

im = array(Image.open('empire.jpg').convert('L'),'f')
print im.shape, im.dtype

控制台输出结果如下所示:

(800, 569, 3) uint8
(800, 569) float32

每行的第一个元组表示图像数组的大小(行、列、颜色通道),紧接着的字符串表示数组元素的数据类型。因为图像通常被编码成无符号八位整数(uint8),所以在第一种情况下,载入图像并将其转换到数组中,数组的数据类型为“uint8”。在第二种情况下,对图像进行灰度化处理,并且在创建数组时使用额外的参数“f”;该参数将数据类型转换为浮点型。关于更多数据类型选项,可以参考图书 [24]。注意,由于灰度图像没有颜色信息,所以在形状元组中,它只有两个数值。

数组中的元素可以使用下标访问。位于坐标 ij,以及颜色通道 k 的像素值可以像下面这样访问:

value = im[i,j,k]

多个数组元素可以使用数组切片方式访问。切片方式返回的是以指定间隔下标访问该数组的元素值。下面是有关灰度图像的一些例子:

im[i,:] = im[j,:]      # 将第 j 行的数值赋值给第 i 行
im[:,i] = 100          # 将第 i 列的所有数值设为100
im[:100,:50].sum()     # 计算前100 行、前 50 列所有数值的和
im[50:100,50:100]      # 50~100 行,50~100 列(不包括第 100 行和第 100 列)
im[i].mean()           # 第 i 行所有数值的平均值
im[:,-1]               # 最后一列
im[-2,:] (or im[-2])   # 倒数第二行

注意,示例仅仅使用一个下标访问数组。如果仅使用一个下标,则该下标为行下标。注意,在最后几个例子中,负数切片表示从最后一个元素逆向计数。我们将会频繁地使用切片技术访问像素值,这也是一个很重要的思想。

我们有很多操作和方法来处理数组对象。本书将在使用到的地方逐一介绍。你可以查阅在线文档或者开源图书 [24] 获取更多信息。

1.3.2 灰度变换

将图像读入 NumPy 数组对象后,我们可以对它们执行任意数学操作。一个简单的例子就是图像的灰度变换。考虑任意函数 f,它将 0...255 区间(或者 0...1 区间)映射到自身(意思是说,输出区间的范围和输入区间的范围相同)。下面是关于灰度变换的一些例子:

from PIL import Image
from numpy import *

im = array(Image.open('empire.jpg').convert('L'))

im2 = 255 - im # 对图像进行反相处理

im3 = (100.0/255) * im + 100 # 将图像像素值变换到100...200 区间

im4 = 255.0 * (im/255.0)**2 # 对图像像素值求平方后得到的图像

第一个例子将灰度图像进行反相处理;第二个例子将图像的像素值变换到 100...200 区间;第三个例子对图像使用二次函数变换,使较暗的像素值变得更小。图 1-4 为所使用的变换函数图像。图 1-5 是输出的图像结果。你可以使用下面的命令查看图像中的最小和最大像素值:

print int(im.min()), int(im.max())

图 1-4:灰度变换示例。三个例子中所使用函数的图像,其中虚线表示恒等变换

图 1-5:灰度变换。对图像应用图 1-4 中的函数:f(x)=255-x 对图像进行反相处理(左);f(x)=(100/255)x+100 对图像进行变换(中);f(x)=255(x/255)2 对图像做二次变换(右)

如果试着对上面例子查看最小值和最大值,可以得到下面的输出结果:

2 255
0 253
100 200
0 255

array() 变换的相反操作可以使用 PIL 的 fromarray() 函数完成:

pil_im = Image.fromarray(im)

如果你通过一些操作将“uint8”数据类型转换为其他数据类型,比如之前例子中的 im3 或者 im4,那么在创建 PIL 图像之前,需要将数据类型转换回来:

pil_im = Image.fromarray(uint8(im))

如果你并不十分确定输入数据的类型,安全起见,应该先转换回来。注意,NumPy 总是将数组数据类型转换成能够表示数据的“最低”数据类型。对浮点数做乘积或除法操作会使整数类型的数组变成浮点类型。

1.3.3 图像缩放

NumPy 的数组对象是我们处理图像和数据的主要工具。想要对图像进行缩放处理没有现成简单的方法。我们可以使用之前 PIL 对图像对象转换的操作,写一个简单的用于图像缩放的函数。把下面的函数添加到 imtool.py 文件里:

def imresize(im,sz):
  """ 使用PIL 对象重新定义图像数组的大小"""
  pil_im = Image.fromarray(uint8(im))

  return array(pil_im.resize(sz))

我们将会在接下来的内容中使用这个函数。

1.3.4 直方图均衡化

图像灰度变换中一个非常有用的例子就是直方图均衡化。直方图均衡化是指将一幅图像的灰度直方图变平,使变换后的图像中每个灰度值的分布概率都相同。在对图像做进一步处理之前,直方图均衡化通常是对图像灰度值进行归一化的一个非常好的方法,并且可以增强图像的对比度。

在这种情况下,直方图均衡化的变换函数是图像中像素值的累积分布函数(cumulative distribution function,简写为 cdf,将像素值的范围映射到目标范围的归一化操作)。

下面的函数是直方图均衡化的具体实现。将这个函数添加到 imtool.py 里:

def histeq(im,nbr_bins=256):
  """ 对一幅灰度图像进行直方图均衡化"""

  # 计算图像的直方图
  imhist,bins = histogram(im.flatten(),nbr_bins,normed=True)
  cdf = imhist.cumsum() # 累积分布函数
  cdf = 255 * cdf / cdf[-1] # 归一化

  # 使用累积分布函数的线性插值,计算新的像素值
  im2 = interp(im.flatten(),bins[:-1],cdf)

  return im2.reshape(im.shape), cdf

该函数有两个输入参数,一个是灰度图像,一个是直方图中使用小区间的数目。函数返回直方图均衡化后的图像,以及用来做像素值映射的累积分布函数。注意,函数中使用到累积分布函数的最后一个元素(下标为 -1),目的是将其归一化到 0...1 范围。你可以像下面这样使用该函数:

from PIL import Image
from numpy import *

im = array(Image.open('AquaTermi_lowcontrast.jpg').convert('L'))
im2,cdf = imtools.histeq(im)

图 1-6 和图 1-7 为上面直方图均衡化例子的结果。上面一行显示的分别是直方图均衡化之前和之后的灰度直方图,以及累积概率分布函数映射图像。可以看到,直方图均衡化后图像的对比度增强了,原先图像灰色区域的细节变得清晰。

图 1-6:直方图均衡化示例。左侧为原始图像和直方图,中间图为灰度变换函数,右侧为直方图均衡化后的图像和相应直方图

图 1-7:直方图均衡化示例。左侧为原始图像和直方图,中间图为灰度变换函数,右侧为直方图均衡化后的图像和相应直方图

1.3.5 图像平均

图像平均操作是减少图像噪声的一种简单方式,通常用于艺术特效。我们可以简单地从图像列表中计算出一幅平均图像。假设所有的图像具有相同的大小,我们可以将这些图像简单地相加,然后除以图像的数目,来计算平均图像。下面的函数可以用于计算平均图像,将其添加到 imtool.py 文件里:

def compute_average(imlist):
  """ 计算图像列表的平均图像"""

  # 打开第一幅图像,将其存储在浮点型数组中
  averageim = array(Image.open(imlist[0]), 'f')

  for imname in imlist[1:]:
    try:
      averageim += array(Image.open(imname))
    except:
      print imname + '...skipped'
  averageim /= len(imlist)

  # 返回uint8 类型的平均图像
  return array(averageim, 'uint8')

该函数包括一些基本的异常处理技巧,可以自动跳过不能打开的图像。我们还可以使用 mean() 函数计算平均图像。mean() 函数需要将所有的图像堆积到一个数组中;也就是说,如果有很多图像,该处理方式需要占用很多内存。我们将会在下一节中使用该函数。

1.3.6 图像的主成分分析(PCA)

PCA(Principal Component Analysis,主成分分析)是一个非常有用的降维技巧。它可以在使用尽可能少维数的前提下,尽量多地保持训练数据的信息,在此意义上是一个最佳技巧。即使是一幅 100×100 像素的小灰度图像,也有 10 000 维,可以看成 10 000 维空间中的一个点。一兆像素的图像具有百万维。由于图像具有很高的维数,在许多计算机视觉应用中,我们经常使用降维操作。PCA 产生的投影矩阵可以被视为将原始坐标变换到现有的坐标系,坐标系中的各个坐标按照重要性递减排列。

为了对图像数据进行 PCA 变换,图像需要转换成一维向量表示。我们可以使用 NumPy 类库中的flatten() 方法进行变换。

将变平的图像堆积起来,我们可以得到一个矩阵,矩阵的一行表示一幅图像。在计算主方向之前,所有的行图像按照平均图像进行了中心化。我们通常使用 SVD(Singular Value Decomposition,奇异值分解)方法来计算主成分;但当矩阵的维数很大时,SVD 的计算非常慢,所以此时通常不使用 SVD 分解。下面就是 PCA 操作的代码:

from PIL import Image
from numpy import *

def pca(X):
  """ 主成分分析:
    输入:矩阵X ,其中该矩阵中存储训练数据,每一行为一条训练数据
    返回:投影矩阵(按照维度的重要性排序)、方差和均值"""

  # 获取维数
  num_data,dim = X.shape

  # 数据中心化
  mean_X = X.mean(axis=0)
  X = X - mean_X

if dim>num_data:
  # PCA- 使用紧致技巧
  M = dot(X,X.T) # 协方差矩阵
  e,EV = linalg.eigh(M) # 特征值和特征向量
  tmp = dot(X.T,EV).T # 这就是紧致技巧
  V = tmp[::-1] # 由于最后的特征向量是我们所需要的,所以需要将其逆转
  S = sqrt(e)[::-1] # 由于特征值是按照递增顺序排列的,所以需要将其逆转
  for i in range(V.shape[1]):
    V[:,i] /= S
else:
  # PCA- 使用SVD 方法
  U,S,V = linalg.svd(X)
  V = V[:num_data] # 仅仅返回前nun_data 维的数据才合理

# 返回投影矩阵、方差和均值
return V,S,mean_X

该函数首先通过减去每一维的均值将数据中心化,然后计算协方差矩阵对应最大特征值的特征向量,此时可以使用简明的技巧或者 SVD 分解。这里我们使用了 range() 函数,该函数的输入参数为一个整数 n,函数返回整数 0...(n-1) 的一个列表。你也可以使用 arange() 函数来返回一个数组,或者使用 xrange() 函数返回一个产生器(可能会提升速度)。我们在本书中贯穿使用range() 函数。

如果数据个数小于向量的维数,我们不用 SVD 分解,而是计算维数更小的协方差矩阵 XXT 的特征向量。通过仅计算对应前 kk 是降维后的维数)最大特征值的特征向量,可以使上面的 PCA 操作更快。由于篇幅所限,有兴趣的读者可以自行探索。矩阵 V 的每行向量都是正交的,并且包含了训练数据方差依次减少的坐标方向。

我们接下来对字体图像进行 PCA 变换。fontimages.zip 文件包含采用不同字体的字符 a 的缩略图。所有的 2359 种字体可以免费下载 2。假定这些图像的名称保存在列表 imlist 中,跟之前的代码一起保存传在 pca.py 文件中,我们可以使用下面的脚本计算图像的主成分:

2免费字体图像库由 Martin Solli 收集并上传(http://webstaff.itn.liu.se/~marso/)。

from PIL import Image
from numpy import *
from pylab import *
import pca

im = array(Image.open(imlist[0])) # 打开一幅图像,获取其大小
m,n = im.shape[0:2] # 获取图像的大小
imnbr = len(imlist) # 获取图像的数目

# 创建矩阵,保存所有压平后的图像数据
immatrix = array([array(Image.open(im)).flatten()
               for im in imlist],'f')

# 执行 PCA 操作
V,S,immean = pca.pca(immatrix)

# 显示一些图像(均值图像和前 7 个模式)
figure()
gray()
subplot(2,4,1)
imshow(immean.reshape(m,n))
for i in range(7):
  subplot(2,4,i+2)
  imshow(V[i].reshape(m,n))

show()

注意,图像需要从一维表示重新转换成二维图像;可以使用 reshape() 函数。如图 1-8 所示,运行该例子会在一个绘图窗口中显示 8 个图像。这里我们使用了 PyLab 库的 subplot() 函数在一个窗口中放置多个图像。

图 1-8:平均图像(左上)和前 7 个模式(具有最大方差的方向模式)

1.3.7 使用pickle模块

如果想要保存一些结果或者数据以方便后续使用,Python 中的 pickle 模块非常有用。pickle模块可以接受几乎所有的 Python 对象,并且将其转换成字符串表示,该过程叫做封装(pickling)。从字符串表示中重构该对象,称为拆封(unpickling)。这些字符串表示可以方便地存储和传输。

我们来看一个例子。假设想要保存上一节字体图像的平均图像和主成分,可以这样来完成:

# 保存均值和主成分数据
f = open('font_pca_modes.pkl', 'wb')
pickle.dump(immean,f)
pickle.dump(V,f)
f.close()

在上述例子中,许多对象可以保存到同一个文件中。pickle 模块中有很多不同的协议可以生成 .pkl 文件;如果不确定的话,最好以二进制文件的形式读取和写入。在其他 Python 会话中载入数据,只需要如下使用 load() 方法:

# 载入均值和主成分数据
f = open('font_pca_modes.pkl', 'rb')
immean = pickle.load(f)
V = pickle.load(f)
f.close()

注意,载入对象的顺序必须和先前保存的一样。Python 中有个用 C 语言写的优化版本,叫做cpickle 模块,该模块和标准 pickle 模块完全兼容。关于 pickle 模块的更多内容,参见pickle 模块文档页 http://docs.python.org/library/pickle.html

在本书接下来的章节中,我们将使用 with 语句处理文件的读写操作。这是 Python 2.5 引入的思想,可以自动打开和关闭文件(即使在文件打开时发生错误)。下面的例子使用 with() 来实现保存和载入操作:

# 打开文件并保存
with open('font_pca_modes.pkl', 'wb') as f:
  pickle.dump(immean,f)
  pickle.dump(V,f)

# 打开文件并载入
with open('font_pca_modes.pkl', 'rb') as f:
  immean = pickle.load(f)
  V = pickle.load(f)

上面的例子乍看起来可能很奇怪,但 with() 确实是个很有用的思想。如果你不喜欢它,可以使用之前的 open 和 close 函数。

作为 pickle 的一种替代方式,NumPy 具有读写文本文件的简单函数。如果数据中不包含复杂的数据结构,比如在一幅图像上点击的点列表,NumPy 的读写函数会很有用。保存一个数组 x 到文件中,可以使用:

savetxt('test.txt',x,'%i')

最后一个参数表示应该使用整数格式。类似地,读取可以使用:

x = loadtxt('test.txt')

你可以从在线文档http://docs.scipy.org/doc/numpy/reference/generated/numpy.loadtxt.html 了解更多内容。

最后,NumPy 有专门用于保存和载入数组的函数。你可以在上面的在线文档里查看关于 save()和 load() 的更多内容。

1.4 SciPy

SciPyhttp://scipy.org/) 是建立在 NumPy 基础上,用于数值运算的开源工具包。SciPy 提供很多高效的操作,可以实现数值积分、优化、统计、信号处理,以及对我们来说最重要的图像处理功能。接下来,本节会介绍 SciPy 中大量有用的模块。SciPy 是个开源工具包,可以从http://scipy.org/Download 下载。

1.4.1 图像模糊

图像的高斯模糊是非常经典的图像卷积例子。本质上,图像模糊就是将(灰度)图像 I 和一个高斯核进行卷积操作:

Iσ = I*Gσ

其中 * 表示卷积操作;Gσ 是标准差为 σ 的二维高斯核,定义为 :

G_\sigma=\frac{1}{2\pi\sigma^2}e^{-(x^2+y^2)/2\sigma^2}

高斯模糊通常是其他图像处理操作的一部分,比如图像插值操作、兴趣点计算以及很多其他应用。

SciPy 有用来做滤波操作的 scipy.ndimage.filters 模块。该模块使用快速一维分离的方式来计算卷积。你可以像下面这样来使用它:

from PIL import Image
from numpy import *
from scipy.ndimage import filters

im = array(Image.open('empire.jpg').convert('L'))
im2 = filters.gaussian_filter(im,5)

上面 guassian_filter() 函数的最后一个参数表示标准差。

图 1-9 显示了随着 σ 的增加,一幅图像被模糊的程度。σ 越大,处理后的图像细节丢失越多。如果打算模糊一幅彩色图像,只需简单地对每一个颜色通道进行高斯模糊:

im = array(Image.open('empire.jpg'))
im2 = zeros(im.shape)
for i in range(3):
  im2[:,:,i] = filters.gaussian_filter(im[:,:,i],5)
im2 = uint8(im2)

在上面的脚本中,最后并不总是需要将图像转换成 uint8 格式,这里只是将像素值用八位来表示。我们也可以使用:

im2 = array(im2,'uint8')

来完成转换。

关于该模块更多的内容以及不同参数的选择,请查看http://docs.scipy.org/doc/scipy/reference/ndimage.html 上 SciPy 文档中的 scipy.ndimage部分。

图 1-9:使用 scipy.ndimage.filters 模块进行高斯模糊:(a)原始灰度图像;(b)使用 σ=2 的高斯滤波器;(c)使用 σ=5 的高斯滤波器;(d)使用 σ=10 的高斯滤波器

1.4.2 图像导数

整本书中可以看到,在很多应用中图像强度的变化情况是非常重要的信息。强度的变化可以用灰度图像 I(对于彩色图像,通常对每个颜色通道分别计算导数)的 x 和 y 方向导数 Ix 和 Iy 进行描述。

图像的梯度向量为∇I = [IxIy]T。梯度有两个重要的属性,一是梯度的大小

\left|\boldsymbol{\nabla I}\right|=\sqrt{{\boldsymbol{I}_x}^2+{\boldsymbol{I}_y}^2}

它描述了图像强度变化的强弱,一是梯度的角度

α=arctan2(IyIx)

描述了图像中在每个点(像素)上强度变化最大的方向。NumPy 中的 arctan2() 函数返回弧度表示的有符号角度,角度的变化区间为 -π...π。

我们可以用离散近似的方式来计算图像的导数。图像导数大多数可以通过卷积简单地实现:

Ix=I*Dx 和 Iy=I*Dy

对于 Dx 和 Dy,通常选择 Prewitt 滤波器:

D_x=\begin{vmatrix}-1&0&1\\-1&0&1\\ -1&0&1\end{vmatrix} 和 D_y=\begin{vmatrix}-1&-1&-1\\0&0&0\\ 1&1&1\end{vmatrix}

或者 Sobel 滤波器:

D_x=\begin{vmatrix}-1&0&1\\-2&0&2\\ -1&0&1\end{vmatrix} 和 D_y=\begin{vmatrix}-1&-2&-1\\0&0&0\\ 1&2&1\end{vmatrix}

这些导数滤波器可以使用 scipy.ndimage.filters 模块的标准卷积操作来简单地实现,例如:

from PIL import Image
from numpy import *
from scipy.ndimage import filters

im = array(Image.open('empire.jpg').convert('L'))

# Sobel 导数滤波器
imx = zeros(im.shape)
filters.sobel(im,1,imx)

imy = zeros(im.shape)
filters.sobel(im,0,imy)

magnitude = sqrt(imx**2+imy**2)

上面的脚本使用 Sobel 滤波器来计算 x 和 y 的方向导数,以及梯度大小。sobel() 函数的第二个参数表示选择 x 或者 y 方向导数,第三个参数保存输出的变量。图 1-10 显示了用 Sobel 滤波器计算出的导数图像。在两个导数图像中,正导数显示为亮的像素,负导数显示为暗的像素。灰色区域表示导数的值接近于零。

图 1-10:使用 Sobel 导数滤波器计算导数图像:(a)原始灰度图像;(b)x 导数图像;(c)y导数图像;(d)梯度大小图像

上述计算图像导数的方法有一些缺陷:在该方法中,滤波器的尺度需要随着图像分辨率的变化而变化。为了在图像噪声方面更稳健,以及在任意尺度上计算导数,我们可以使用高斯导数滤波器:

Ix=I*Gσx 和 Iy=I*Gσy

其中,Gσx和 Gσy 表示 Gσ 在 x 和 y 方向上的导数,Gσ 为标准差为 σ 的高斯函数。

我们之前用于模糊的 filters.gaussian_filter() 函数可以接受额外的参数,用来计算高斯导数。可以简单地按照下面的方式来处理:

sigma = 5 # 标准差

imx = zeros(im.shape)
filters.gaussian_filter(im, (sigma,sigma), (0,1), imx)

imy = zeros(im.shape)
filters.gaussian_filter(im, (sigma,sigma), (1,0), imy)

该函数的第三个参数指定对每个方向计算哪种类型的导数,第二个参数为使用的标准差。你可以查看相应文档了解详情。图 1-11 显示了不同尺度下的导数图像和梯度大小。你可以和图 1-9 中做相同尺度模糊的图像做比较。

图 1-11:使用高斯导数计算图像导数:x 导数图像(上),y 导数图像(中),以及梯度大小图像(下);(a)为原始灰度图像,(b)为使用 σ=2 的高斯导数滤波器处理后的图像,(c)为使 用 σ=5 的高斯导数滤波器处理后的图像,(d)为使用 σ=10 的高斯导数滤波器处理后的图像

1.4.3 形态学:对象计数

形态学(或数学形态学)是度量和分析基本形状的图像处理方法的基本框架与集合。形态学通常用于处理二值图像,但是也能够用于灰度图像。二值图像是指图像的每个像素只能取两个值,通常是 0 和 1。二值图像通常是,在计算物体的数目,或者度量其大小时,对一幅图像进行阈值化后的结果。你可以从 http://en.wikipedia.org/wiki/Mathematical_morphology 大体了解形态学及其处理图像的方式。

scipy.ndimage 中的 morphology 模块可以实现形态学操作。你可以使用 scipy.ndimage 中的measurements 模块来实现二值图像的计数和度量功能。下面通过一个简单的例子介绍如何使用它们。

考虑在图 1-12a3 里的二值图像,计算该图像中的对象个数可以通过下面的脚本实现:

3这个图像实际上是图像“分割”后的结果。如果你想知道该图像是如何创建的,可以查看 9.3 节。

from scipy.ndimage import measurements,morphology

# 载入图像,然后使用阈值化操作,以保证处理的图像为二值图像
im = array(Image.open('houses.png').convert('L'))
im = 1*(im<128)

labels, nbr_objects = measurements.label(im)
print "Number of objects:", nbr_objects

上面的脚本首先载入该图像,通过阈值化方式来确保该图像是二值图像。通过和 1 相乘,脚本将布尔数组转换成二进制表示。然后,我们使用 label() 函数寻找单个的物体,并且按照它们属于哪个对象将整数标签给像素赋值。图 1-12b 是 labels 数组的图像。图像的灰度值表示对象的标签。可以看到,在一些对象之间有一些小的连接。进行二进制开(binary open)操作,我们可以将其移除:

# 形态学开操作更好地分离各个对象
im_open = morphology.binary_opening(im,ones((9,5)),iterations=2)

labels_open, nbr_objects_open = measurements.label(im_open)
print "Number of objects:", nbr_objects_open

binary_opening() 函数的第二个参数指定一个数组结构元素。该数组表示以一个像素为中心时,使用哪些相邻像素。在这种情况下,我们在 y 方向上使用 9 个像素(上面 4 个像素、像素本身、下面 4 个像素),在 x 方向上使用 5 个像素。你可以指定任意数组为结构元素,数组中的非零元素决定使用哪些相邻像素。参数 iterations 决定执行该操作的次数。你可以尝试使用不同的迭代次数 iterations 值,看一下对象的数目如何变化。你可以在图 1-12c 与图 1-12d 中查看经过开操作后的图像,以及相应的标签图像。正如你想象的一样,binary_closing() 函数实现相反的操作。我们将该函数和在 morphology 和 measurements 模块中的其他函数的用法留作练习。你可以从 scipy.ndimage 模块文档 http://docs.scipy.org/doc/scipy/reference/ndimage.html 中了解关于这些函数的更多知识。

图 1-12:形态学示例。使用二值开操作将对象分开,然后计算物体的数目:(a)为原始二值图像;(b)为对应原始图像的标签图像,其中灰度值表示物体的标签;(c)为使用开操作后的二值图像;(d)为开操作后图像的标签图像

1.4.4 一些有用的SciPy模块

SciPy 中包含一些用于输入和输出的实用模块。下面介绍其中两个模块:io 和 misc

  1. 读写.mat文件

    如果你有一些数据,或者在网上下载到一些有趣的数据集,这些数据以 Matlab 的 .mat 文件格式存储,那么可以使用 scipy.io 模块进行读取。

    data = scipy.io.loadmat('test.mat')
    

    上面代码中,data 对象包含一个字典,字典中的键对应于保存在原始 .mat 文件中的变量名。由于这些变量是数组格式的,因此可以很方便地保存到 .mat 文件中。你仅需创建一个字典(其中要包含你想要保存的所有变量),然后使用 savemat() 函数:

    data = {}
    data['x'] = x
    scipy.io.savemat('test.mat',data)
    

    因为上面的脚本保存的是数组 x,所以当读入到 Matlab 中时,变量的名字仍为 x。关于scipy.io 模块的更多内容,请参见在线文档http://docs.scipy.org/doc/scipy/reference/io.html

  2. 以图像形式保存数组

    因为我们需要对图像进行操作,并且需要使用数组对象来做运算,所以将数组直接保存为图像文件 4 非常有用。本书中的很多图像都是这样的创建的。

    imsave() 函数可以从 scipy.misc 模块中载入。要将数组 im 保存到文件中,可以使用下面的命令:

    from scipy.misc import imsave
    imsave('test.jpg',im)
    

    scipy.misc 模块同样包含了著名的 Lena 测试图像:

    lena = scipy.misc.lena()
    

    该脚本返回一个 512×512 的灰度图像数组。

4所有 Pylab 图均可保存为多种图像格式,方法是点击图像窗口中的“保存”按钮。

1.5 高级示例:图像去噪

我们通过一个非常实用的例子——图像的去噪——来结束本章。图像去噪是在去除图像噪声的同时,尽可能地保留图像细节和结构的处理技术。我们这里使用 ROF(Rudin-Osher-Fatemi)去噪模型。该模型最早出现在文献 [28] 中。图像去噪对于很多应用来说都非常重要;这些应用范围很广,小到让你的假期照片看起来更漂亮,大到提高卫星图像的质量。ROF 模型具有很好的性质:使处理后的图像更平滑,同时保持图像边缘和结构信息。

ROF 模型的数学基础和处理技巧非常高深,不在本书讲述范围之内。在讲述如何基于 Chambolle 提出的算法 [5] 实现 ROF 求解器之前,本书首先简要介绍一下 ROF 模型。

一幅(灰度)图像 I 的全变差(Total Variation,TV)定义为梯度范数之和。在连续表示的情况下,全变差表示为:

J(\boldsymbol{I})=\int\left|\nabla\boldsymbol{I}\right|\text{dx}            (1.1)

在离散表示的情况下,全变差表示为:

J(\boldsymbol{I})=\sum_{\text{x}}\left|\nabla\boldsymbol{I}\right|

其中,上面的式子是在所有图像坐标 x=[x, y] 上取和。

在 Chambolle 提出的 ROF 模型里,目标函数为寻找降噪后的图像 U,使下式最小:

\min_U\left|\left|\boldsymbol{I}-\boldsymbol{U}\right|\right|^2+2\lambda J(\boldsymbol{U}),

其中范数 ||I-U|| 是去噪后图像 U 和原始图像 I 差异的度量。也就是说,本质上该模型使去噪后的图像像素值“平坦”变化,但是在图像区域的边缘上,允许去噪后的图像像素值“跳跃”变化。

按照论文 [5] 中的算法,我们可以按照下面的代码实现 ROF 模型去噪:

from numpy import *

def denoise(im,U_init,tolerance=0.1,tau=0.125,tv_weight=100):
  """ 使用A. Chambolle(2005)在公式(11)中的计算步骤实现Rudin-Osher-Fatemi(ROF)去噪模型

    输入:含有噪声的输入图像(灰度图像)、U 的初始值、TV 正则项权值、步长、停业条件

    输出:去噪和去除纹理后的图像、纹理残留"""

  m,n = im.shape # 噪声图像的大小

  # 初始化
  U = U_init
  Px = im # 对偶域的x 分量
  Py = im # 对偶域的y 分量
  error = 1

while (error > tolerance):
  Uold = U

  # 原始变量的梯度
  GradUx = roll(U,-1,axis=1)-U # 变量U 梯度的x 分量
  GradUy = roll(U,-1,axis=0)-U # 变量U 梯度的y 分量

  # 更新对偶变量
  PxNew = Px + (tau/tv_weight)*GradUx
  PyNew = Py + (tau/tv_weight)*GradUy
  NormNew = maximum(1,sqrt(PxNew**2+PyNew**2))

  Px = PxNew/NormNew # 更新x 分量(对偶)
  Py = PyNew/NormNew # 更新y 分量(对偶)

  # 更新原始变量
  RxPx = roll(Px,1,axis=1) # 对x 分量进行向右x 轴平移
  RyPy = roll(Py,1,axis=0) # 对y 分量进行向右y 轴平移

  DivP = (Px-RxPx)+(Py-RyPy) # 对偶域的散度
  U = im + tv_weight*DivP # 更新原始变量

  # 更新误差
  error = linalg.norm(U-Uold)/sqrt(n*m);

  return U,im-U # 去噪后的图像和纹理残余

在这个例子中,我们使用了 roll() 函数。顾名思义,在一个坐标轴上,它循环“滚动”数组中的元素值。该函数可以非常方便地计算邻域元素的差异,比如这里的导数。我们还使用了linalg.norm() 函数,该函数可以衡量两个数组间(这个例子中是指图像矩阵 U和 Uold)的差异。我们将这个 denoise() 函数保存到 rof.py 文件中。

下面使用一个合成的噪声图像示例来说明如何使用该函数:

from numpy import *
from numpy import random
from scipy.ndimage import filters
import rof

# 使用噪声创建合成图像
im = zeros((500,500))
im[100:400,100:400] = 128
im[200:300,200:300] = 255
im = im + 30*random.standard_normal((500,500))

U,T = rof.denoise(im,im)
G = filters.gaussian_filter(im,10)

# 保存生成结果
from scipy.misc import imsave
imsave('synth_rof.pdf',U)
imsave('synth_gaussian.pdf',G)

原始图像和图像的去噪结果如图 1-13 所示。正如你所看到的,ROF 算法去噪后的图像很好地保留了图像的边缘信息。

图 1-13:使用 ROF 模型对合成图像去噪:(a)为原始噪声图像;(b)为经过高斯模糊的图像(σ=10);(c)为经过 ROF 模型去噪后的图像

下面看一下在实际图像中使用 ROF 模型去噪的效果:

from PIL import Image
from pylab import *
import rof

im = array(Image.open('empire.jpg').convert('L'))
U,T = rof.denoise(im,im)

figure()
gray()
imshow(U)
axis('equal')
axis('off')
show()

经过 ROF 去噪后的图像如图 1-14c 所示。为了方便比较,该图中同样显示了模糊后的图像。可以看到,ROF 去噪后的图像保留了边缘和图像的结构信息,同时模糊了“噪声”。

图 1-14:使用 ROF 模型对灰度图像去噪:(a)为原始噪声图像;(b)为经过高斯模糊的图像(σ=5);(c)为经过 ROF 模型去噪后的图像

练习

  1. 如图 1-9 所示,将一幅图像进行高斯模糊处理。随着 σ 的增加,绘制出图像轮廓。在你绘制出的图中,图像的轮廓有何变化?你能解释为什么会发生这些变化吗?

  2. 通过将图像模糊化,然后从原始图像中减去模糊图像,来实现反锐化图像掩模操作(http://en.wikipedia.org/wiki/Unsharp_masking)。反锐化图像掩模操作可以实现图像锐化效果。试在彩色和灰度图像上使用反锐化图像掩模操作,观察该操作的效果。

  3. 除了直方图均衡化,商图像是另一种图像归一化的方法。商图像可以通过除以模糊后的图像I/(IGσ) 获得。尝试使用该方法,并使用一些样本图像进行验证。

  4. 使用图像梯度,编写一个在图像中获得简单物体(例如,白色背景中的正方形)轮廓的函数。

  5. 使用梯度方向和大小检测图像中的线段。估计线段的长度以及线段的参数,并在原始图像中重新绘制该线段。

  6. 使用 label() 函数处理二值化图像,并使用直方图和标签图像绘制图像中物体的大小分布。

  7. 使用形态学操作处理阈值化图像。在发现一些参数能够产生好的结果后,使用 morphology 模块里面的 center_of_mass() 函数寻找每个物体的中心坐标,将其在图像中绘制出来。

代码示例约定

从第 2 章起,我们假定 PIL、NumPy 和 Matplotlib 都包括在你所创建的每个文件和每个代码例子的开头:

from PIL import Image
from numpy import *
from pylab import *

这种约定使得示例代码更清晰,同时也便于读者理解。除此之外,我们使用 SciPy 模块时,将会在代码示例中显式声明。

一些纯化论者会反对这种将全体模块导入的方式,坚持如下使用方式:

import numpy as np
import matplotlib.pyplot as plt

这种方式能够保持命名空间(知道每个函数从哪儿来)。因为我们不需要 PyLab 中的 NumPy 部分,所以该例子只从 Matplotlib 中导入 pyplot 部分。纯化论者和经验丰富的程序员们知道这些区别,他们能够选择自己喜欢的方式。但是,为了使本书的内容和例子更容易被读者接受,我们不打算这样做。

请读者注意。