• 基于图像处理的汽车牌照的识别作者:陈秋菊指导老师:李方洲 (温州师范学院 物理与电子信息学院 325027) 摘要:以一幅汽车牌照的识别为例,具体介绍了车牌自动识别的原理。整个处理过程分为预处理、边缘提取、车牌...

               

     

    基于图像处理的汽车牌照的识别

    作者:陈秋菊

    指导老师:李方洲

                 (温州师范学院 物理与电子信息学院  325027)

     

    摘要:以一幅汽车牌照的识别为例,具体介绍了车牌自动识别的原理。整个处理过程分为预处理、边缘提取、车牌定位、字符分割、字符识别五大模块,用MATLAB软件编程来实现每一个部分,最后识别出汽车牌照。在研究的同时对其中出现的问题进行了具体分析,处理。寻找出对于具体的汽车牌照识别过程的最好的方法。

     

    关键词汽车牌照 车牌提取 字符分割  字符识别

     

    The vehicle license recognition based on the image processing

             Author:Chen Qiuju

                Tutor:Li Fangzhou

    (School of Physics and Electronic Information   Wen Zhou Normal College  325027)

    Abstract With one vehicle license recognition, the principle of the automobile License recognition is introduced .This process was divided into pre-process, edge extraction, vehicle license location, character division and character recognition, which is implemented separated by using MATLAB. The license is recognized at last. At the same time, the problems are also analyzed

    And solved in the process. The best method of recognition to the very vehicle license is found.

     

    Keywords: vehicle license   vehicle license location   character segmentation

    Character  recognition

                 

     

     

     

     

     

     

     

     

    1. 引言

    1.1 选题意义

       汽车牌照自动识别系统是以汽车牌照为特定目标的专用计算机视觉系统,是计算机视觉和模式识别技术在智能交通领域应用的重要研究课题之一,是实现交通管理智能化的重要环节,它可广泛应用于交通流量检测,交通控制与诱导,机场、港口、小区的车辆管理,不停车自动收费,闯红灯等违章车辆监控以及车辆安全防盗等领域,具有广阔的应用前景。目前,发达国家LPR汽车牌照识别技术License Plate Recognition LPR,简称“车牌通”)系统在实际交通系统中已成功应用,而我国的开发应用进展缓慢,车牌识别系统基本上还停留在实验室阶段基于这种现状还有它广阔的应用前景,目前对汽车车牌的识别研究就有了深远的意义。

                课题组成

    汽车车牌的识别过程主要包括车牌定位、字符车牌分割和车牌字符识别三个关键环节。其识别流程如下:

    字符识别  

    字符分割 

    车牌定位  

    边缘提取  

    图像预处理  

    原始图像   

                                                                      

     


    原始图像  :由数码相机或其它扫描装置拍摄到的图像

    图像预处理:对动态采集到的图像进行滤波,边界增强等处理以克服图像干扰

    边缘提取  :通过微分运算,2值化处理,得到图像的边缘

    车牌定位  :计算边缘图像的投影面积,寻找峰谷点,大致确定车牌位置,再计算此连通域内的宽高比,剔除不在域值范围内的连通域。最后得到的便为车牌区域。

    字符分割  :利用投影检测的字符定位分割方法得到单个的字符

    字符识别  :利用模板匹配的方法与数据库中的字符进行匹配从而确认出字符,得到最后的汽车牌照,包括英文字母和数字。

    本文以一幅汽车图像为例,结合图像处理各方面的知识,利用MATLAB编程,实现了从

    车牌的预处理到字符识别的完整过程。各部分的处理情况如下:

     

     

     

     

     

     

     

     

     

     

     

     

    2.预处理及边缘提取

     

                              1 汽车原图

       图像在形成、传输或变换过程中,受多种因素的影响,如:光学系统失真、系统噪声、暴光不足或过量、相对运动等,往往会与原始景物之间或图像与原始图像之间产生了某种差异,这种差异称为降质或退化。这种降质或退化对我们的处理往往会造成影响。因此在图像处理之前必须进行预处理,包括去除噪音,边界增强,增加亮度等等。

    因为噪声主要是一些含高频的突变成分,因此可以通过一个低通滤波器来消除图像中包含的噪声,并使低频成分得到增强。滤波的方式有两种,一种是空间域滤波,一种是频率域滤波。在空间域,常见的滤波方式有两种方式,均值滤波和中值滤波。空间域滤波主要有巴特沃斯滤波器。在车牌边缘提取之前,两种滤波方式均采用了。并与未进行滤波的边缘进行比较。以下是经处理后的一些图片。

            2  经均值滤波后提取的边缘图像

           3 经巴特沃斯低通滤波后提取的边缘图像

              4 未滤波直接提取出的边缘信息

        

              5 经高通滤波器增强后得到的边缘图像

       对比以上几幅图片,图2的边缘太粗,而图3的边缘已经模糊掉了。图5中包含的噪声太多,图4未经滤波直接提取出的边缘图像最清晰,所包含的有用信息最多。分析这种情况产生的原因,归纳起来主要有以下方面:

      1、原始图像清晰度比较高,从而简化了预处理

      2、图像的平滑处理会使图像的边缘信息受到损失,图像变得模糊

      3、图像的锐化可以增强图像中物体的边缘轮廓,但同时也使一些噪声得到了增强

      综上所述,结合MATLAB实验过程,得出不是每一种图像处理之初都适合滤波和边界增强。本次汽车车牌的识别,为了保存更多的有用信息,经过多次比较,选择图4作为后期处理的依据。边缘的提取采用的是梯度算子,因为其实现过程比较简单,所以在此不多加赘述了。

      提取出的边缘含有多个灰度值,要进行二值化处理,选择一个合适的域值。经多次比较,选取域值T=70,对于灰度值大于T的赋值为255,小于T的赋值为0。经过处理后的图像如下所示:

     

                6 二值化后的边缘图像

    结合后期分割得到的车牌图7,二值化后的图像在后期的识别中并不会提高车牌的识别率,因此不采用二值化的图像来进行识别,因此后面的处理依然使用图4

    3.车牌提取

        经过边缘提取得到的图像,车牌区域在水平方向灰度面积值具有明显频繁的跳变,在垂直方向上的面积投影则出现峰--峰的特性。根据这种峰谷特点,自动检测车牌位置峰点检测的车牌区域定位方法并对初步定位后的车牌进一步使用微定位技术该方法包括三部分: (1) 车牌的横向定位 (2)车牌的纵向定位 (3) 车牌的微定位。

       汽车本身具有一定的特点一般情况下牌照都挂在缓冲器上或附近处于车牌照图像的下半部分本次分割的主要意图是缩小牌照搜索范围大致确定出牌照的位置。对如图4所示的汽车边缘图像f ( x y) 我们首先进行水平方向一阶差分运算

    g ( i j) = | f ( i j) - f ( i j + 1) |

    其中i = 1 2 3 xw - 1 j = 0 1 2 3 yw - 1 ,其中xwyw分别为图像的行数和列数。然后对水平差分图像的像素沿水平方向累加产生一个投影表T( i) 如图7所示。

                   7 汽车边缘图像的水平面积投影图

    一般对应于车牌位置的投影值T( i) 较大而在车牌上下行附近的投影值较小均有谷点存在。只要能找到这两个谷点就能大致确定出车牌照的位置缩小车牌搜索范围。由图4可以看出,车牌下方的横栏处的T(i)值应该是最大的,而车牌位置就在其附近。根据这些特定,可定出车牌位置大概在320~350行之间。

      类似的方法得出汽车边缘图像的垂直面积投影图

                  8 汽车边缘图像的垂直面积投影图

       同上可初步得到汽车牌照的列位置在120~210之间。大致确定的牌照位置如下图。

                     

                      9 粗略定位出的汽车牌照

    对初步确定出来的牌照进行微定位,所谓微定位法就是对基本定位后的车牌图像进

    行局部分析以进一步确定字符范围缩减车牌的左、右和上、下边界这有利于后续的牌照字符处理。具体实现如下: (1) 由于车牌近似为一个矩形上下边缘近似为一条直线通过简单的灰度变化分析就可以再次定位车牌图像的上下边界这种情况适合于倾斜度较小的车牌对于倾斜程度较大的牌照来说在其横向定位之前就应该利用相关的技术进行车牌的矫正(例如Hough 变换技术) (2)确定左边界: 从左向右扫描 ,当遇到灰度值大于设定值60之后,停止扫描。上边界也是利用这种方式得到。这样就得到首字符的起始位置。再利用牌照的大小,宽高比一般都是固定的这些先验知识,就可以确定出牌照的具体位置。本设计中采用的车牌,其宽高比为1:3。从而确定出汽车牌照的具体位置。最后提取出的汽车牌照如下图:

    10 二值化的汽车牌照

                     

    11 未进行二值化的汽车牌照

    4.字符分割

    在汽车牌照自动识别过程中,字符分割有承前启后的作用。它在前期牌照定位的基础

    上进行字符的分割,然后再利用分割的结果进行字符识别。字符识别的算法很多,常采用垂直面积投影法来实现。面积投影法的公式如下:

                

    由于字符块在竖直方向上的投影必然在字符间或字符内的间隙处取得局部最小值,并且这个位置应满足车牌的字符书写格式、字符尺寸和其他一些条件的限制。下图是图10在垂直方向上的面积投影图。从图形中我们很直观的看出投影值中出现了8条间隙, 6个字母中间的间隙只有5个,还有三个间隙是字符间的。有字符的列其灰度值比较高,无字符的则相对比较低。依据这一点,再结合图10的特征,很容易得到每个字符的起始终止位置。第一个字符:1-10 第二个字符:10-18 第三个字符:28-41 第四个字符:42-48第五个字符:60-68 第六个字符:68-78

      

                         12 车牌垂直方向上的面积投影图

      将图10按照上面的分析行数不变,列数分为六组,分别影射到六个不同的数组中。又因为在字符的模式识别中,其模板大小统一,因此得到的六个数组必须变换其大小,均统一成26×14的形式。分割出来的六个字符如下所示,分别命名为M1.jpgM2.jpgM3.jpgM4.jpgM5.jpgM6jpg并用imwrite函数写入图像文件夹中,以便在后期处理中可以直接进行调用。

                                    

                           12 分割出来的六个字符图像

    一般分割出来的字符要进行进一步的处理,以满足下一步字符识别的需要。因为图像中含有许多燥声,这在预处理的图像中已经看出来了。因此必须进行滤波,然后归一化,二值处理。使其最后得到的图像与标准模板一样。只含有两种灰度值,黑与白。但是对于车牌的识别,并不需要这么多的处理就已经可以达到正确识别的目的。在此简化了处理过程,未经滤波归一化,直接进行后期处理。

    5.字符识别

    字符的识别目前用于车牌字符识别(OCR)中的算法主要有基于模板匹配的OCR算法以及基于人工神经网络的OCR算法。基于模板匹配的OCR的基本过程是:首先对待识别字符进行二值化并将其尺寸大小缩放为字符数据库中模板的大小,然后与所有的模板进行匹配,最后选最佳匹配作为结果。用人工神经网络进行字符识别主要有两种方法:一种方法是先对待识别字符进行特征提取,然后用所获得的特征来训练神经网络分类器。识别效果与字符特征的提取有关,而字符特征提取往往比较耗时。因此,字符特征的提取就成为研究的关键。另一种方法则充分利用神经网络的特点,直接把待处理图像输入网络,由网络自动实现特征提取直至识别。模板匹配的主要特点是实现简单,当字符较规整时对字符图像的缺损、污迹干扰适应力强且识别率相当高。综合模板匹配的这些优点我们将其用为车牌字符识别的主要方法。

     

    字符识别的算法如下:

                       

    读入字符图像

       

     

     


    与模板库中的字母逐一进行相关运算

                

     

     

    寻找相关度最大值所对应的模板

     

     


                     

     

    输出此模板所对应的值

     

     

     

     


      因此在字符识别之前必须把模板库设置好。汽车拍照的字符一般有七个,大部分车牌第一位是汉字,通常代表车辆所属省份,或是军种、警别等有特定含义的字符简称;紧接其后的为字母与数字。车牌字符识别与一般文字识别在于它的字符数有限,汉字共约50多个,大写英文字母26个,数字10个。所以建立字符模板库也极为方便。本次设计所识别的车牌只有字母与数字。为了实验方便,结合本次设计所选汽车牌照的特点,只建立了3个字母与3个数字的模板。其他模板设计的方法与此相同。

    5.1模板设计

    分析字符分割得到的图像以及其他车牌图像中字符的特点,将模板大小定为26×14

    背景为黑色,代表灰度值 0,字符边缘为白色。代表灰度值255

    设计过程如下:

    1.  用画图工具先画出PFM103等几个字符的图像。并分别保存为m11.jpgm12.jpg

    M1.jpgm14.jpgm15.jpgm16.jpg。根据画图的经验其大小应略大于26×14,以利于后面的处理。所得到的字符均为黑字白底。字体为方正姚体,大小16号。

    所画出的图形如下:

    2.     1中获得的粗略图像用MATLAB进行处理。处理的方面包括a.去除边框,利用二维数组的映射关系实现,去掉开始行()与结束行()即所有外围像素点b.提取模板的边缘。其实现方法与前面的汽车边缘提取的方法相同,也是利用梯度算子。C.尺寸大小变换,切割为标准的26×14格式。

    3.     整个处理过程结束后,再用imwrite 函数写入图像库中,作为标准模板使用。

    经过处理后的标准模板形式如下:

                                 

    模板的设计过程,共进行了两次。作为一幅jpg或者bmp形式的图片,其中包含了许许

    多多的像素点,各象素点的值也不一样,到底如何来确定这些像素点的值呢?模板设置无疑成为一个难点。于是便投机取巧,将本身分割得到的字符保存起来直接作为模板使用。这种方法明显不具有科学性。图像识别中模板匹配法是将未知的东西与已知的东西来进行比较。因此模板库中的字符应该是已知的,是预先设置好的是不可更改的,而不是在实验过程中来得到。第二次便是用以上的画图工具来一步一步逼近实现模板设置。实验结果显示,这种方法简单易行。

    5.2识别过程

      字符识别中模板匹配方法是实现离散输入模式分类的有效途径之一,其实质是度量输入与样本之间的某种相似性,取相似性最大者为输入模式所属类别,它根据字符的直观形象抽取特征,用相关匹配原理进行识别,即是将输入字符与标准字符在一个分类器中进行匹配。车牌字符相关匹配算法如下:

       输入字符用输入函数X表示,标准模板用函数T表示,它们的大小均为26×14。将未知的模式逐个与模板匹配,求出其相似度。

      

              

          

      其中,M=26N=14Ti  指第i个模板。X×Ti 是指矩阵对应元素相乘。如果把XT都归一化为1”0”,则上式表示标准模板与待识图像上对应点均为1”的数目与标准模板上1”点的数目之比。

      如果 maxSi>λ 则判定XTi,否则拒识,这里λ为拒识域值。本次实验中直接取相关最大值为判定域值。

      MATLAB来实现字符的识别。其中有两个很重要的函数可以直接调用:corr2 函数和max函数。

      Corr2 函数直接用来计算图像拒阵AB的相关系数rAB的大小及数据类型必须相同,得到的结果r是双精度标量。其调用方法如下:

                    r=corr2(AB)

      因此本设计中的字符识别只需调用此函数,即将分割出来的字符与设置好的模板一一进行相关运算,然后寻找出它们中的最大相关值。max函数就是用来选择几个数值的大小值。其调用方式如下:

                    [bc]=max(a(:))

      其中b返回的是比较后得到的最大值,c是最大值所对应的元素位置。

    同时还自编了一个识别函数result.m用来返回识别结果。程序代码如下:

      function c=result(H)

    M1=imread('M1.jpg')

    M2=imread('M2.jpg')

    M3=imread('M3.jpg')

    M4=imread('M4.jpg')

    M5=imread('M5.jpg')

    M6=imread('M6.jpg')

    M1=double(M1)

    M2=double(M2)

    M3=double(M3)

    M4=double(M4)

    M5=double(M5)

    M6=double(M6)

    d=zeros(6)

    d(1)=corr2(HM1)

    d(2)=corr2(HM2)

    d(3)=corr2(HM3)

    d(4)=corr2(HM4)

    d(5)=corr2(HM5)

    d(6)=corr2(HM6)

    [De]=max(d(:))

    switch e

        case 1

            c='P'

        case 2

            c='F'

        case 3

            c='M'

        case 4

            c=1

        case 5

            c=0

        case 6

            c=3

        otherwise

    end

       M1 M2 M3 M4 M5 M6为标准模板对应的二维数组。函数在开始运行时,从主函数里得到待识别的字符,此字符与模板进行相关运算,同时找出相关运算后的最大值及对应的函数位置。然后对最大位置进行判别,返回所对应的字符给主程序。并在MATLAB运行窗口显示出来。

      编写result.m函数大大简化了计算量,每个字符在识别的时候直接调用此函数,避免了重新编程浪费的时间及空间。这是本次车牌识别中的一大亮点

      运行结果,MATLAB的运行窗口显示出的车牌号码为:PFM103。完成了准确识别车牌的目的。

    6.总结:

      6.1 设计过程说明

      在汽车车牌识别的整个过程中,查找了很多资料,综合了各方面的信息。车牌实现的每一步都有许多的方法,各种方法都有其优劣,但是对于具体的图像处理,并不是每一种理论在实践中都可以实现,即使实现了也很难说哪一种方法最合适,还得在具体的实验中比较选择。以车牌预处理为例,经过滤波和边界增强处理后提取出的车牌效果明显没有未处理的效果理想。第二点在程序调试的过程中要耐心的检查每一个错误。测试结果表明,本设计有以下几条优点:

    1.     充分利用MATLAB中已有的函数库。整个程序设计简单易行

    2.     识别准确率高

    62设计工具说明:

      车牌识别程序设计能够得以顺利完成。在很大程度上得利于MATLAB这套软件, MATLAB功能强大,它包括数值计算和符号计算,并且计算结果和编程可视化。这为编程调试创造了一个便利的环境。作为图像处理最适用的工具之一,其突出的特点是它包含一个图像处理工具包,这个工具包由一系列支持图像处理操作的函数组成。所支持的图像处理操作有:图像的几何操作、邻域和区域操作、图像变换、图像恢复与增强、线性滤波和滤波器设计、变换(DCT 变换等) 、图像分析和统计、二值图像操作等。在图像的显示方面,MATLAB 提供了图像文件读入函数imread ( ) 用来读取如: bmp tif jpgpcxtiff jpeghdf xwd 等格式图像文件图像写出函数imwrite () 还有图像显示函数image ( ) imshow() 等等。这些,都使编程效率大为提高。

    7.讨论

       本设计已成功达到了车牌识别的目的,而且准确率也很高。但是在整个设计过程中,发现了几个值得参考的算法,也试图用这种算法来实现车牌识别,但种种原因,而未采用。现在次简单讨论一下:

      6.1.8个方向的模板的边缘提取算法

    本设计在边缘提取之时,使用的是梯度算子。这种算法的优点是简单易行,但是因为它不包含边缘的方向,因此对噪声不够敏感。目前,在图像处理方面使用得最多的是一种可抗噪声的Sobel算法。它定义了8个方向的模板。

                   

    通常物体的边缘是连续且光滑的,而噪声是随机的。 在任一边缘点附近沿边缘的走向总能找到另一边缘点,且这两边缘点之间的灰度差及方向差都不可能很大。但是噪声则不同,一般情况下,沿任一噪声点的方向(通过上述模板运算得到)不太可能找到与其灰度差及方向差都很小的噪声点。正是利用这一基本思想,本算法能将实际的边缘点与噪声点区分开来。

    加权领域平均算法来进行滤波处理

     62 加权领域平均算法来进行滤波处理

     由实验我们可以看出,一般的滤波器在对图像进行噪声滤除的同时对图像中的细节部分有不同程度的破坏,都不能达到理想的效果。但是采用加权的邻域平均算法对图像进行噪声滤除不仅能够有效地平滑噪声还能够锐化模糊图像的边缘。 加权的邻域平均算法的基本思想是: 在一个邻域内除了可以利用灰度均值外灰度的上偏差和下偏差也能够提供某些局部信息。算法的计算公式描述如下f (x y ) 表示原始图像g (x y ) 为平滑后点(x y ) 的灰度值V x y 表示以点(x y ) 为中心的邻域该邻域包含N 个象素m (x y ) 表示邻域V x y 内的灰度均值。NI表示邻域内大于平均值的像素个数,Ng表示小于平均值的像素个数,而N0表示等于平均值的像素个数。则修正的邻域平均法由下式给出:

     

    m - Aı m lN l > max{N g N 0}

    g(xy)=   m + Aı m gN g > max{N l N 0}            (1)

    m    else

     

     (1)(1) A为修正系数取值范围为01其大小反映V x y 中的边缘状况。 以上是我认为在图像处理中比较有价值的两点,有兴趣的可以上网查阅相关的资料。

    【参考文献】

    [1] 霍宏涛.数字图像处理.机械工业出版社,2003.5

    [2] 陈桂明、张明照、戚红雨.应用MATLAB语言处理数字信号与数字图像。科学出版社,2000

    [3] 郎锐.数字图象处理学Visual C++实现.北京希望电子出版社,2002.12

    [4] 刘露、强.汽车牌照自动识别技术初探,中国公路网,2003-09-26

    [5] 周妮娜、王敏、黄心汉、吕雪峰、万国红.车牌字符识别的预处理算法.计算机工程与应用,2003(15)

    [6] 佘新平、朱 立.一种具有抗噪声干扰的图像边缘提取算法的研究2001-6-4

    [7] 苑玮琦、伞晓钟. 一种汽车牌照多层次分割定位方法,2004 Vol.9 No.4 P.239-243

    [8] 许志影李晋平.MATLAB极其在图像处理中的应用.计算机与现代化,2004(4)

    [9] 董慧颖、曹仁帅.汽车牌照自动识别系统中字符分割算法研究.沈阳工业学院学报2003(12)Vol.22 No.4

    [10] 崔 江王友仁.车牌自动识别方法中的关键技术研究.计算机测量与控制,2003.11(4)

    [11] 王年、李婕、任彬、汪炳权.多层次汽车车牌照定位分割方法. 安徽大学学报,1999(6)Vol.23.No.2

    [12] 马俊莉莫玉龙王明祥.一种基于改进模板匹配的车牌字符识别方法.小型微型计算机系统Vol.23.No.2

      

     

     

     

     

     

     

     

     
    展开全文
  • 形态学图像处理中,输入的是图像,输出的是从图像中提取出来的属性,分割是该方向上的另一步骤。 分割将图像细分为构成它的子区域或物体。细分的程度取决于要解决的问题。当感兴趣的物体或区域已经被检测出来了,那...

    参考:https://blog.csdn.net/mary_0830/article/details/89597672
    https://blog.csdn.net/Dujing2019/article/details/90203492


    形态学图像处理中,输入的是图像,输出的是从图像中提取出来的属性,分割是该方向上的另一步骤。
    分割将图像细分为构成它的子区域或物体。细分的程度取决于要解决的问题。当感兴趣的物体或区域已经被检测出来了,那就停止分割。
    单色图像的分割算法通常基于灰度值的两个基本性质:不连续性和相似性。第一种是以灰度变换(即灰度的突变)为基础分割图像,如边缘检测;第二种是根据事先预定义的相似性准则把图像分成相似的区域,如阈值分割、区域生长、区域分裂和聚合。

    图像分割的数学描述
    令R表示一幅图像占据的整个空间区域。将图像分割视为把R分为n个子区域R1,R2,…,Rn的过程,满足:
    在这里插入图片描述
    (a):分割必须是完全的,即每个像素都必须在一个区域内;
    (b):一个区域中的点以某种预定义的方式来连接(这些点必须是4连接或8连接的);
    ©:各个区域必须不相交;
    (d):分割后的区域中的像素必须满足的属性(如,若Ri中的所有像素都有相同的灰度级,则Q(Ri) = TRUE),其中属性可以是灰度,颜色,纹理等;
    (e):两个邻接区域Ri和Rj在属性Q的意义上必须是不同的。

    一、点、线和边缘检测

    本节主要介绍以灰度局部剧烈变化(灰度不连续性)检测为基础的分割方法。
    对一幅图中灰度的突变,局部变化可以用微分来检测。因为变化短促,可以用一阶微分和二阶微分描述。
    数字函数的导数可用差分来定义。一维函数f(x)在点x处一阶导数的近似:
    在这里插入图片描述
    一维函数f(x)在点x处二阶导数的近似:
    在这里插入图片描述
    对二维图像函数f(x,y),为表示一致性,使用偏微分,将处理沿两个空间轴的偏微分

    有关一阶导数和二阶导数的结论:
    (1) 一阶导数通常在图像中产生较粗的边缘
    (2) 二阶导数对精细细节,如细线、孤立点和噪声有较强的响应;
    (3) 二阶导数在灰度斜坡和灰度台阶过渡处会产生双边缘响应
    (4) 二阶导数的符号可用于确定边缘的过渡是从亮到暗(负二阶导数)还是从暗到亮(正二阶导数)

    计算图像中每个像素位置处的一阶导数和二阶导数的另一种方法是使用空间滤波器
    对3*3滤波器掩膜来说,导数是计算模板系数与被该模板覆盖的区域中的灰度值的乘积之和。即模板在该区域中心点处的响应如下:
    在这里插入图片描述
    zk为像素的灰度,wk为对应位置的系数。
    其中,3 x 3的空间滤波器为:
    在这里插入图片描述

    1.1 孤立点的检测

    点的检测应以二阶导数为基础。
    孤立点常嵌在常数区域(或图像中亮度基本不变的区域)中,孤立点的灰度将完全不同与其周围像素。
    点检测通常使用拉普拉斯算子进行检测。第三章介绍的四种拉普拉斯模板如下:
    在这里插入图片描述
    原理:如果在某个点处,该模板的响应R的绝对值超过了某个指定阈值,则说明在模板中心位置(x,y)处的该点已被检测到。输出图像中,这样的点标注为1,其他所有的点标注为0,从而产生一幅二值图像。公式描述如下:
    在这里插入图片描述
    其中,g是输出图像,T为一个非负阈值。R为上述介绍过的,即模板的系数与其覆盖区域的灰度值的乘积之和,也叫作模板的响应。
    当模板的中心位于一个孤立点时,模板的响应最强,而在亮度不变的区域响应为0。

    实现
    若T已经给出

    g = abs(imfilter(double(f),w))>=T;  
    f:输入图像;
    w:点检测模板;
    T:指定的非负阈值;
    g:包含检测点的输出图像
    

    示例:图像中的孤立点检测

    f = imread('test_pattern_with_single_pixel.tif');
    w = [1 1 1; 1 -8 1; 1 1 1];
    g = abs(imfilter(double(f),w));
    T = max(g(:)); % 将滤波后的图像g的最大值设定为阈值T,很明显,g中不可能有比T值更大的点。
    g = g >= T;
    
    figure;
    subplot(1,2,1); imshow(f); title('原图');
    subplot(1,2,2); imshow(g); title('点检测图像');
    

    在这里插入图片描述
    由图可以看出,原图方框框的位置,有一个几乎看不见的小黑点,通过点检测的方法,在检测后的图像中,就能很明显的看到该点的存在。
    点检测的另一种方法是在大小为m*n的所有邻域内寻找一些点,这些点的最大值和最小值之差超过了T。可以使用ordfilt2函数(排序)来实现。

    g = imsubtract(ordfilt2(f,m*n,ones(m,n)),ordfilt2(f,1,ones(m,n)));
    g = g >= T;
    

    1.2 线检测

    线检测会更复杂一些。对于线检测,二阶导数将产生更强的响应,并产生比一阶导数更细的线。
    ==注:==必须适当处理二阶导数的双线效应
    线检测模板如下:
    在这里插入图片描述
    第一个模板对水平方向的线有最佳响应;第二个模板对45°方向的线有最佳响应;第三个模板对垂直方向的线有最佳响应;第四个模板对**-45°方向**的线有最佳响应。每个模板的首选方向用一个比其他方向更大的系数(如2)加权。每个模板中系数之和为0,表明恒定灰度区域中的响应为0。
    :有关正负45°模板中,不同的书有不同的说法,均与所选的坐标系及其方向有关,只要自己能搞清楚到底想要的是图像哪个方向的边缘即可,根据自己想要提取的对角的方向选择合适的模板即可,而不必纠结到底哪个模板是+45°,哪个是-45°。

    若对检测图像中由给定模板定义的方向上的所有线感兴趣,则只须简单对图像运用这个模板,并对结果的绝对值进行阈值处理。留下的点是有最强响应的点,这些点与模板的方向最为接近,且组成了只有一个像素宽的线。

    示例:寻找图像中所有宽度为1像素、方向为-45°的线。

    f = imread('wirebond_mask.tif');
    w = [2 -1 -1;-1 2 -1;-1 -1 2];
    g = imfilter(double(f),w);
    
    gtop = g(1:120,1:120); % 获取左上角的图像
    gtop = pixeldup(gtop,4); % 将图像放大4*4倍
    
    gbot = g(end -119:end,end -119:end);  % 获取右下角的图像
    gbot = pixeldup(gbot,4);  % 将图像放大4*4倍
    
    g1 = abs(g); 
    T = max(g1(:));  % 将阈值设为g1的最大值
    g2 = g1 >= T; % 对模板滤波的结果进行阈值处理
    
    figure;
    subplot(2,3,1),imshow(f),title('(a)原图');
    subplot(2,3,2),imshow(g,[]),title('(b)45°方向处理后的图像');
    subplot(2,3,3),imshow(gtop,[]),title('(c)左上角放大的图像');
    subplot(2,3,4),imshow(gbot,[]),title('(d)右下角放大的图像');
    subplot(2,3,5),imshow(g1,[]),title('(e)滤波后绝对值的图像');
    subplot(2,3,6),imshow(g2),title('(f)阈值处理后的图像');
    

    在这里插入图片描述
    由于一个figure中显示的图太多,故有些图中的细节可能显示不完整,若想要观察每个图的具体细节,单独显示要观察的图就好。
    从图b中可以看到,图中有明显的黑线和白线,这就是二阶导数的双线效应,中等灰度表示0,较暗的灰色调表示负值,较量的色调表示正值。放大图像会看的更明显。对图b求绝对值,可以消除图中的负值,如图e所示。
    从图c和图d中的直线段相比,图d中的更亮,这是因为图a中右下角的线段宽度为1,而左上方线段的宽度不是。

    1.3 边缘检测

    边缘检测是基于灰度突变来分割图像的方法

    1.3.1 边缘模型

    边缘模型根据它们的灰度剖面分为:台阶模型、斜坡模型和屋顶模型。
    在这里插入图片描述
    台阶边缘:理想型,1像素的距离上出现两个灰度级间的理想过渡。
    斜坡边缘:实际中,数字图像都存在被模糊且带有噪声的边缘,因此更接近斜坡边缘,斜坡的斜度与边缘的模糊程度成正比。一个边缘点是斜坡中包含的任何点,一条边缘线段是一组已连接起来的这样的点。
    屋顶边缘:屋顶模型是通过一个区域的线的模型,屋顶边缘的基底(宽度)由该线的宽度和尖锐度决定。
    在这里插入图片描述
    结论:
    (1) 一阶导数的幅度可用于检测图像中的某个点处是否存在一个边缘;
    (2) 二阶导数的符号可用于确定一个边缘像素是位于该边缘的暗侧还是位于该边缘的亮侧,为“正”在暗侧,为“负”在亮侧;
    围绕一条边缘的二阶导数的两个附加性质:
    a、对图像中的每条边缘,二阶导数生成两个值(一个不希望的特点);
    b、二阶导数的零交叉点(零灰度轴和二阶导数极值间的连线的交点)可用于定位粗边缘的中心

    边缘检测的基本步骤:
    在这里插入图片描述
    (1) 平滑滤波:由于梯度计算易受噪声影响,因此第一步是用滤波去噪。但是,降低噪声的平滑能力越强,边界强度的损失越大;
    (2) 锐化滤波:为了检测边界,必须确定某点邻域中灰度的变换。锐化操作加强了存在有意义的灰度局部变化位置的像素点;
    (3) 边缘判定:在图像中存在许多梯度不为零的点,但对于特定应用,不是所有点都有意义。这就要求我们根据具体情况选择和去除处理点,具体方法包括二值化处理和过零检测等;
    (4) 边缘连接:将间断的边缘连接成有意义的完整边缘,同时去除假边缘。主要方法是Hough变换

    1.3.2 基本边缘检测

    图像梯度
    图像f的(x,y)处边缘的强度和方向用梯度表征,梯度用∇f表示,用向量来定义:
    在这里插入图片描述
    该向量指出了f在位置(x,y)处的最大变化率的方向
    ∇f的大小(长度/赋值)表示为M(x,y)。
    在这里插入图片描述
    平方和平方根需要大量的计算开销,故常用的是用绝对值来近似梯度的大小(幅值)。
    在这里插入图片描述
    它是梯度向量方向变化率的值,M(x,y)常称为梯度图像:gx,gy和M(x,y)都是与原图像大小相同的图像。
    梯度向量的方向由对于x轴度量的角度给出:
    在这里插入图片描述
    任意点(x,y)处一个边缘的方向与该点处梯度向量的方向α(x,y)正交。

    梯度算子(一阶导数)
    常见的梯度算子的模板如下:
    在这里插入图片描述
    Roberts交叉梯度算子
    具有对角优势的二维模板之一。
    在这里插入图片描述
    Prewitt算子
    检测水平和竖直方向。
    在这里插入图片描述
    Sobel算子
    与Prewitt模板类似,检测水平和竖直方向,但中心系数上使用一个权值2:
    在这里插入图片描述
    中心位置处使用2可以平滑图像。

    用于检测对角线方向边缘的Prewitt和Sobel算子
    在这里插入图片描述
    小结:
    1、Roberts算子边缘定位精确度高,但易丢失一部分边缘,由于图像没经过平滑处理,故不具备抑制噪声的能力。该算子对具有陡峭边缘且含噪声少的图像效果较好
    2、Sobel算子和Prewitt算子都考虑了邻域信息,相当于对图像先做加权平滑处理,不同的是平滑部分的权值有差异,故对噪声具有一定的抑制能力,检测出的边缘容易出现多像素宽度。Prewitt模板比Sobel模板实现更简单,但Sobel模板抑制噪声能力更强。
    3、当边缘检测的目的是突出主要边缘并尽可能保持连接时,实践中通常对图像先进行平滑处理,在边缘检测完之后,对边缘结果再进行阈值处理

    1.3.3 更先进的边缘检测技术

    这些技术是通过考虑如图像噪声和边缘本身特性等因素来改善简单的边缘检测方法。
    高斯拉普拉斯算子(LoG)
    拉普拉斯算子是一个二阶导数,对噪声具有无法接受的敏感性,而且其幅值会产生双线效应,另外,边缘方向的不可检测也是其的缺点之一,故一般不以其原始形式作用域边缘检测。
    G(x,y)为高斯函数为:
    在这里插入图片描述
    其中σ是标准差。用高斯函数卷积模糊一幅图像,图像的模糊程度是由σ决定的。
    高斯拉普拉斯算子为∇2G(x,y):
    在这里插入图片描述
    常见的LoG模板
    在这里插入图片描述
    实际中使用的是该模板的负值
    算法步骤:
    1、对输入图像使用n*n高斯低通滤波器进行滤波;
    2、计算滤波后图像的拉普拉斯;
    在这里插入图片描述
    3、寻找g(x,y)的零交叉来确定f(x,y)中边缘的位置。
    高斯低通滤波器大小的选择,即n值应是大于等于6σ的最小奇整数。

    坎尼边缘检测器(Canny算子)
    前面介绍的几种基于微分方法的边缘检测算法,都是只有在图像不含噪声或者先通过平滑去噪的前提下才能正常使用。
    在图像边缘检测中,抑制噪声和边缘精确定位是无法同时满足的。Canny算子力图在抗噪声干扰和精确定位之间寻求最佳这种方案。
    Canny方法基于三个基本目标:
    1、低错误率。所有边缘都应该被找到,并且应该没有伪响应。即检测到的边缘必须尽可能是真实的边缘。数学上,就是使信噪比SNR尽量大,输出信噪比越大,错误率越少。
    2、边缘点应被很好地定位。已定位边缘必须尽可能接近真实边缘。即由检测器标记为边缘的点和真实边缘的中次年之间的距离应该最小。
    3、单一的边缘点响应。对于真实的边缘点,检测器仅应返回一个点。即真实边缘周围的局部最大数应该是最小的。这意味着在仅存一个单一边缘点的位置,检测器不应指出多个边缘像素。

    算法步骤:
    1、用一个高斯滤波器(n*n大小)平滑输入图像;
    2、用一阶偏导的有限差分来计算梯度幅值图像和角度图像;
    3、对梯度幅值图像应用非最大抑制;
    4、用双阈值处理和连接分析来检测并连接边缘

    1.3.4 边缘检测方法的matlab实现

    使用edge函数可以方便地实现上述几种边缘检测方法。该函数的作用是检测灰度图像中的边缘,并返回一个带有边缘信息的二值图像,其中黑色表示背景,白色表示原图像的边缘部分。
    参考:edge函数
    默认情况下,edge函数使用Sobel算子进行边缘检测。
    基于梯度算子的边缘检测

    BW = edge(I,type,thresh,direction,'nothinning')
    [BW,threshOut] = edge(I,type,thresh,direction,'nothinning')
    

    I:需要检测边缘的输入图像,为灰度图像或二值图像。
    type:表示梯度算子的种类。
    在这里插入图片描述
    thresh:敏感度阈值参数,任何灰度值低于此阈值的边缘将不会被检测到。其默认值是空矩阵[],此时算法会自动计算阈值。
    direction:指定了我们感兴趣的边缘方向,edge函数将只检测到direction中指定方向的边缘,默认的是both。
    :如果选择 Roberts 算子,则 ‘horizontal’ 方向实际上检测到与水平方向成 135° 角的边缘,‘vertical’ 方向检测到与水平方向成 45° 角的边缘。
    在这里插入图片描述
    可选参数’nothinning’,指定时可以通过跳过边缘细化算法来加快算法运行的速度。默认时候,这个参数是’thinning’,即进行边缘细化。
    BW:返回的二值图像,其中0(黑色)为背景,1(白色)为边缘部分。
    threshOut:返回阈值。若指定阈值了,则threshOut = thresh,若不指定阈值,则返回自动确定的阈值。

    基于高斯-拉普拉斯算子的边缘检测

    BW = edge(I,'log',thresh,sigma)
    

    ‘log’:表示高斯-拉普拉斯算子。
    thresh:敏感度阈值参数。如果将thresh设为0,则输出的边缘图像将包含围绕所有物体的闭合的轮廓线,因为这样的运算会包括输入图像中所有的过零点。
    sigma:指定生成高斯滤波器所使用的标准差。默认为2。滤镜大小n x n,n的计算方法为n = ceil(sigma*3)*2 + 1。
    BW为返回的二值图像,其中0(黑色)为背景,1(白色)为边缘部分。
    零交叉检测器:

    BW = edge(I,'zerocross',thresh,h) % 使用指定的滤波器 h 检测边缘。
    

    基于Canny算子的边缘检测

    BW = edge(I,'canny',thresh,sigma)
    

    Canny算子是edge函数中最强的边缘检测算子
    ‘Canny’:表示canny算子
    thresh:敏感度阈值参数,其默认值是空矩阵[]。Canny算法的敏感度阈值是一个列向量,因为需要为算法指定阈值的上下限。忽略边缘强度低于下阈值的所有边缘,保留边缘强度高于上阈值的所有边缘。在指定阈值矩阵时,第1个元素是阈值下限,第2个元素为阈值上限。如果只指定一个阈值元素,则这个直接指定的值会被作为阈值上限,而它与0.4的乘积会被作为阈值下限。如果阈值参数没有指定,算法会自行确定敏感度阈值的上下限。
    sigma:指定生成平滑使用的高斯滤波器的标准差。默认为1。滤镜大小n x n,n的计算方法为n = ceil(sigma*3)*2 + 1。

    示例:使用Sobel算子来提取边缘

    f = imread('building.tif');
    [gv,t] = edge(f,'sobel','vertical'); % 提取垂直边缘
    % 主要边缘是垂直边缘(倾斜的边缘也有垂直分量和水平分量),故倾斜边缘也能被检测到
    t  % 默认阈值为0.0516
    gv1 = edge(f,'sobel',0.15,'vertical');  %  自行指定一个较高的阈值,将弱一些的边缘去掉,如0.15
    
    gboth = edge(f,'sobel',0.15); % 水平和垂直边缘
    
    wneg45 = [-2 -1 0;-1 0 1;0 1 2];  % 用于提取45°方向的边缘
    w45 = [0 1 2;-1 0 1;-2 -1 0];  % 用于提取-45°方向的边缘
    gneg45 = imfilter(double(f),wneg45,'replicate');
    g45 = imfilter(double(f),w45,'replicate');
    T =0.3*max(abs(gneg45(:)));
    gneg45 = gneg45>=T;
    g45 = g45 >= T;
    
    figure;
    subplot(2,3,1); imshow(f); title('原图');
    subplot(2,3,2); imshow(gv);title('默认阈值垂直的sobel图像');
    subplot(2,3,3); imshow(gv1); title('阈值0.15垂直的sobel图像');
    subplot(2,3,4); imshow(gboth); title('阈值为0.15垂直、水平的sobel图像');
    subplot(2,3,5); imshow(gneg45);title('45°边缘图像');
    subplot(2,3,6);imshow(g45); title('-45°边缘图像');
    

    在这里插入图片描述
    示例:五种边缘检测算法的比较

    f = imread('building.tif');
    bw1 = edge(f,'sobel');
    bw2 = edge(f,'prewitt');
    bw3 = edge(f,'roberts');
    bw4 = edge(f,'log');
    bw5 = edge(f,'canny');
    
    figure;
    subplot(2,3,1); imshow(f); title('原图');
    subplot(2,3,2); imshow(bw1); title('sobel检测');
    subplot(2,3,3); imshow(bw2); title('prewitt检测');
    subplot(2,3,4); imshow(bw3); title('roberts检测');
    subplot(2,3,5); imshow(bw4); title('log检测');
    subplot(2,3,6); imshow(bw5); title('canny检测');
    

    在这里插入图片描述
    分析:
    1、从边缘定位的精度看:
    Roberts算子和Log算子定位精度较高
    Roberts算子简单直观,Log算子利用二阶导数零交叉特性检测边缘。但Log算子只能获得边缘位置信息,不能得到边缘的方向等信息。
    2、从对不同方向边缘的响应看
    从对边缘方向的敏感性而言,Sobel算子、Prewitt算子检测斜向阶跃(灰度突变)边缘效果较好,Robets算子检测水平和垂直边缘效果较好。Log算子不具备边缘方向检测能力。
    Sobel算子可以提供最精确的边缘方向估计
    3、从去噪能力看
    Roberts和Log算子定位精度虽然较高,但受噪声影响大。
    Sobel算子和Prewitt算子模板相对较大因为去噪能力较强,具有平滑作用,能滤除一些噪声,去掉部分伪边缘,但同时也平滑了真正的边缘,这也正是其定位精度不高的原因。

    从总体效果来衡量,Canny算子给出了一种边缘定位精确性和抗噪声干扰性的较好折中。

    当有成本和速度限制时,通常使用阈值梯度边缘检测方法。当追求边缘质量时,通常会选择Log和Canny算子,Canny算子效果更好。

    1.3.5 Hough变换

    理想情况下,边缘检测应该仅产生位于边缘上的像素集合。实际上,由于噪声、不均匀照明引起的边缘间断,以及其他引入灰度值虚假的不连续的影响,这些像素并不能完全描述边缘特性。因此,一般是在边缘检测后紧跟连接算法,将边缘像素组合成有意义的边缘或区域边界。
    Hough变换是一个非常重要的检测间断点边界形状的方法。通过将图像坐标空间变换到参数空间,来实现直线和曲线的拟合。

    原理:
    直角坐标参数空间
    在x-y坐标空间中,经过点(xi,yi)的直线表示为yi = axi + b,a为斜率,b为截距。
    通过点(xi,yi)的直线有无数条,对应的a和b也不尽相同。
    若将xi和yi看作常数,将a和b看作变量,从x-y空间变换到a-b参数空间。则点(xi,yi)处的直线变为b = -xia + yi。x-y空间的另一点(xj,yj)处的直线变为b = -xja + yj
    x-y空间中的点在a-b空间中对应一条直线,若点(xi,yi)和(xi,yi)在x-y空间共线,则在a-b空间对应的两直线相交于一点(a’,b’)。反之,在a-b空间相交于同一点的所有直线,在x-y空间都有共线的点与之对应。
    在这里插入图片描述
    具体计算:
    将a-b空间视为离散的。建立二维累加数组A(a,b),第一维是x-y空间中直线斜率的范围,第二维是直线截距范围。二维累加数组A也常被称为Hough矩阵
    初始化A(a,b)为0。
    对x-y空间的每个前景点(xi,yi),将a-b空间的每个a带入b = -xia + yi,计算对应的b。
    每计算出一对(a,b),对应的A(a,b) = A(a,b) + 1。
    所有计算结束后,在a-b空间找最大的A(a,b),即峰值。峰值所对应的(a’,b’)参考点就是原图像中共线点数目最多的直线方程的参数。
    在这里插入图片描述
    求出参考点后,整个目标的边界就可以确定了。

    极坐标参考空间
    使用直角坐标表示直线时,当直线为一条垂直直线或接近垂直直线时,该直线的斜率为无限大或接近无限大,故在a-b空间中无法表示,因此要在极坐标参考空间解决这一问题。
    直线的法线表示为:xcosθ + ysinθ = ρ
    ρ:直线到原点的垂直距离,取值范围为[-D,D],D为一幅图像中对角间的最大距离;θ:x轴到直线垂直线的角度,取值范围为[-90°,90°]。
    极坐标中的Hough变换,是将图像x-y空间坐标变换到ρ-θ参数空间中。x-y空间中共线的点变换到ρ-θ空间后,都相交于一点。不同于直角坐标的是,x-y空间共线的点(xi,yi)和(xj,yj)映射到ρ-θ空间是两条正弦曲线,相交于点(ρ‘,θ’)。
    具体计算时,也要在ρ-θ空间建立一个二维数组累加器A。除了ρ和θ的取值范围不同,其余与直角坐标类似,最后得到的最大的A所对应的(ρ,θ)。
    在这里插入图片描述
    Hough也能处理其他任意形状的函数。

    matlab实现
    通过Hough变换在二值图像中检测直线的3个步骤:

    • 用hough函数执行霍夫变换,得到霍夫矩阵;
    • 用houghpeaks函数在霍夫矩阵中寻找峰值点;
    • 用houghlines函数在原二值图像上提取线段

    1、霍夫变换——hough
    hough函数对一幅二值图像执行Hough变换,得到Hough矩阵。调用格式为:

     [H,theta,rho] = hough(f)   
     或 [H,theta,rho] = hough(f,'ThetaRes',val1,'RhoRes',val2) 
    

    其中,H是霍夫变换矩阵,theta和rho是ρ和θ值向量。f是二值图像,val1是0-90的标量,指定沿θ轴单位区间的间距(默认为1);val2可取(0,norm(size(f))区间上的实数,指定沿ρ轴单位区间的间距(默认为1)。
    2、寻找峰值——houghpeaks
    在Hough矩阵中寻找指定数目的峰值点。调用格式为:

    peaks = houghpeaks(H,NumPeaks)
    或 peaks = houghpeaks(H,NumPeaks,‘Threshold’,val1,'NHoodSize',val2)
    

    H:从hough函数得到的霍夫矩阵
    Numpeaks:要寻找的峰值数目,默认为1
    Threshold:峰值的阈值,只有大于该阈值的点才被认为是可能的峰值,val1可以从0到Inf变换,默认是0.5 x max(H(:))
    NHoodSize:在每次检测出一个峰值后,NHoodSize指出在该峰值周围需要清零的邻域信息。以向量[M,N]的形式给出,其中M、N均为正的奇数,默认为≥size(H)/50的最小奇数。
    peaks:是一个Q*2的矩阵,每行的两个元素分别为某一峰值点在Hough矩阵中的行、列索引,Q为找到的峰值点的数目。
    3、提取直线段——houghlines
    根据Hough矩阵的峰值检测结果提取直线段。调用格式为:

    lines = houghlines(f,theta,rho,peaks)
    或lines = houghlines(f,theta,rho,peaks,'FillGap',val1,'MinLength',val2)
    

    theta 和 rho 是函数 hough 返回的向量。peaks 是由 houghpeaks 函数返回的矩阵。
    FillGap:线段合并的阈值,如果对应于Hough矩阵某一单元格(相同的ρ和θ)的2个线段之间的距离小于FillGap,则合并为1个直线段,默认值为20。
    MinLength:检测直线段的最小长度阈值,如果检测出的直线段长度大于MinLength,则保留,丢弃所有长度小于MinLength,默认值为40。
    lines:是一个结构体数组,数组长度是找到的直线条数,而每一个数组元素的内部结构如下:

    含义
    point1 直线段的端点1坐标
    point2 直线段的端点2坐标
    theta 对应在霍夫矩阵中的θ
    rho 对应在霍夫矩阵中的ρ

    示例:使用霍夫变换连接边缘

    I  = imread('building.tif'); 
    
    % 获取图像边缘(二值图像)
    BW = edge(I,'canny');    
    
    % 执行Hough变换并显示Hough矩阵
    [H,T,R] = hough(BW);
    imshow(H,[],'XData',T,'YData',R,'InitialMagnification','fit');
    xlabel('\theta'), ylabel('\rho');
    axis on, axis normal, hold on;
    
    % 在Hough矩阵中寻找前5个大于Hough矩阵中最大值0.3倍的峰值
    P  = houghpeaks(H,5,'threshold',ceil(0.3*max(H(:))));
    x = T(P(:,2)); y = R(P(:,1));  % 由行、列索引转换成实际坐标
    plot(x,y,'s','color','white');  % 在Hough矩阵图像中标出峰值位置
    
    % 找到并绘制直线
    lines = houghlines(BW,T,R,P,'FillGap',20,'MinLength',10);  % 合并距离小于20的线段,丢弃所有长度小于10的直线段
    figure, imshow(I), hold on
    max_len = 0;
    for k = 1:length(lines)  % 依次标出各条直线段
       xy = [lines(k).point1; lines(k).point2];  % xy中包含一条线段的起点和终点坐标
       plot(xy(:,1),xy(:,2),'LineWidth',2,'Color','green');  % 线为绿色
    
       % 绘制线段端点
       plot(xy(1,1),xy(1,2),'x','LineWidth',2,'Color','yellow');  % 起点为黄色
       plot(xy(2,1),xy(2,2),'x','LineWidth',2,'Color','red');  % 终点为红色
    
       % 确定最长的线段
       len = norm(lines(k).point1 - lines(k).point2);
       if ( len > max_len)
          max_len = len;
          xy_long = xy;
       end
    end
    
    % 高亮显示最长的线段
    plot(xy_long(:,1),xy_long(:,2),'LineWidth',2,'Color','cyan'); 
    

    在这里插入图片描述
    方框框出来的为峰值在参数空间中的位置
    在这里插入图片描述
    将提取出来的直线段叠加在原图中,并对直线加以标注。

    :Hough只能处理二值图像,一般在执行变换前需要在图像上执行边缘检测。

    三、阈值处理

    图像阈值处理在图像分割应用中处于核心地位

    3.1 基础知识

    一般情况下,一幅图像由前景部分和背景部分构成,我们感兴趣的一般的是前景部分,所以一般用阈值将前景和背景分割开来,使我们感兴趣的图像像素值为1,不感兴趣的为0;有时一张图我们会有几个不同的感兴趣区域(不在同一个灰度区域),这时我们可以用多个阈值进行分割。
    单阈值处理
    在这里插入图片描述
    T为阈值。此时T适用于整个图像,该处理称为全局阈值处理。
    双阈值处理
    在这里插入图片描述
    T1和T2分别为不同的阈值,a、b和c是三个不同的灰度值。
    在这里插入图片描述
    灰度阈值的成功与否直接关系到可区分直方图模式(即不同的感兴趣区域,可能不止一个)的波谷的宽度和深度,波谷越宽且越深,则越易进行阈值处理。
    影响波谷特性的关键因素:1)波峰间的间隔(波峰离得越远,分离这些模式的机会越大);2)图像中的噪声内容(模式的波峰随噪声的增加而展宽);3)物体和背景的相对尺寸;4)光源的均匀性;5)图像反射特性的均匀性。

    3.2 基本的全局阈值处理

    物体和背景像素的灰度分布十分明显时,可用适用于整个图像的单个阈值。下面介绍为每幅图像自动估计阈值的算法。算法描述如下:
    1)为全局阈值T选择一个初始估计值(建议初始估计值设为图像最大灰度值和最小灰度值的中间值);
    2)使用阈值T分割该图像。将产生两组像素:G1由灰度值大于T的所有像素组成,G2由灰度值小于T的所有像素组成。
    3)对G1和G2的像素分别计算平均灰度值m1和m2;
    4)计算新阈值:T = 1/2(m1 + m2);
    5)重复步骤2~4,直到在连续的重复中,T的差异比预先设定的参数小为止。

    使用imbinarize进行阈值分割图像

    BW = imbinarize(I,T)
    
    f = imread('scanned-text-grayscale.tif');
    count = 0;
    T = mean2(f); % 初始化阈值T,求整个矩阵像素的平均值
    done = false;
    while ~done
    count  = count+1;  % 记录迭代次数
    g = f>T;
    Tnext = 0.5*(mean(f(g))+mean(f(~g)));
    done = abs(T-Tnext)<0.5;
    T = Tnext;
    end
    
    g = imbinarize(f,T/255);  % 需要将阈值归一化
    subplot(1,3,1); imshow(f); title('原图');
    subplot(1,3,2); imhist(f); title('原图的直方图');
    subplot(1,3,3); imshow(g); title('全局阈值分割的结果');
    

    在这里插入图片描述

    3.3 用Otsu方法的最佳全局阈值处理

    Otsu是一种基于图像直方图的方法。
    图像的阈值处理也可以理解为把像素分配个两个或多个类。Otsu采用最大化类间方差的方法,来确定图像分割的最佳阈值k。
    在这里插入图片描述
    σB2 (k)为类间方差;P1(k)为阈值k时,像素被分到类C1中的概率;m1(k)为分配到类C1的像素的平均灰度值;mG为整个图像的平均灰度值;P1(k)和m1(k)与1类似。
    当设置的方差越大,则越接近正确分割一幅图像的阈值。上式中的k就是我们所要寻找的最佳阈值,当k不唯一时,则将所有的最佳阈值进行取平均值即可。

    实现
    graythresh 使用 Otsu 方法计算全局图像阈值。调用格式如下

    T = graythresh(I)
    [T,EM] = graythresh(I)  
    % I为灰度图像;T为归一化阈值;
    % EM为阈值的有效性度量,以范围 [0,1] 内的正标量形式返回。下界只能由具有单一灰度级的图像获得,上界只能由二值图像获得。一般来说,高EM值说明灰度分成两类的可能性比较高。
    

    示例

    f = imread('scanned-text-grayscale.tif');
    T = mean2(f); % 初始化阈值T,求整个矩阵像素的平均值
    done = false;
    while ~done
    g = f>T;
    Tnext = 0.5*(mean(f(g))+mean(f(~g)));
    done = abs(T-Tnext)<0.5;
    T =Tnext;
    end
    T
    g = imbinarize(f,T/255);  % 需要将阈值归一化
    
    [T1,EM] = graythresh(f);
    EM
    T1*255
    g1 = imbinarize(f,T1);
    subplot(2,2,[1,2]); imshow(f); title('原图');
    subplot(2,2,3); imshow(g); title('基本全局阈值分割的结果');
    subplot(2,2,4); imshow(g1); title('Ostu最佳阈值处理');
    

    在这里插入图片描述
    从这个小例子中,我们看不出来基本全局阈值处理和Ostu最佳阈值处理的区别,但对于某些图像的直方图没有明显的波谷的情况,Ostu方法的阈值处理效果会更好。

    3.4 用图像平滑改善全局阈值处理

    噪声会对图像的阈值处理产生很大的影响。当噪声不能在源头减少时,在阈值处理之前可以先对图像进行平滑处理,以便后续进行更好的阈值处理。

    f = imread('cygnusloop_Xray_original.tif');
    f = uint8(f);
    T = graythresh(f);
    g = imbinarize(f,T);
    
    w = fspecial('average',5);
    fa = imfilter(f,w,'replicate');
    
    Ta = graythresh(fa);
    ga =imbinarize(fa,Ta);
    figure;
    subplot(2,3,1); imshow(f); title('带有噪声的图像');
    subplot(2,3,2); imhist(f); title('噪声图像直方图');
    subplot(2,3,3); imshow(g); title('Ostu对噪声图像的处理结果');
    subplot(2,3,4); imshow(fa); title('去噪图像');
    subplot(2,3,5); imhist(fa); title('去噪图像直方图');
    subplot(2,3,6); imshow(ga); title('Ostu对去噪图像的处理结果');
    

    在这里插入图片描述
    由原图的直方图可以看出,几乎观察不到要进行分类的不同类别的波谷,因此分割出来的结果就不是很好,对图像使用均值滤波器先进行平滑处理(去噪)后再进行阈值处理,分割的效果就很好了。

    3.5 利用边缘改进全局阈值处理

    对于边界明显的图像,前景和背景面积悬殊,但是整体灰度值相近,不易用otsu直接找出正确的阈值,可以使用边缘改进的阈值处理。
    边缘改进的阈值处理:主要是处理那些位于或接近物体和背景间边缘的像素,使得这些像素分离开的操作。
    算法描述:
    1)计算一副边缘图像(梯度和拉普拉斯都可以)。
    2)指定一个阈值T。通常通过百分比指定阈值,常选择边缘图像直方图中百分比相对高(如90%)的值作为阈值(边缘图像直方图中的较高的值表示边缘比较强烈的部分)。
    3)用步骤2中的阈值对步骤1中产生的图像进行阈值处理,产生一副二值图像gT(x,y)。将选择出来的“强”边缘像素的二值图像作为模板。
    4)使用gT(x,y)与原图像f(x,y)相乘,得到一幅新图像,相当于使原图像中的物体和背景对比的更明显。
    5)对新图像进行阈值分割,如Ostu方法。
    代码是参考网上的,本人没有弄的很清楚
    示例1:基于梯度边缘信息改进全局阈值处理

    f = double(imread('noisy-image.tif'));
    f = imnoise(f,'gaussian');
    subplot(231),imshow(f),title('原图像');
    subplot(232),imhist(f),title('原图像直方图');
    % 计算梯度
    sx = fspecial('sobel');
    sy = sx';
    gx = imfilter(f,sx,'replicate');
    gy = imfilter(f,sy,'replicate');
    grad = sqrt(gx.*gx+gy.*gy);
    grad = grad/max(grad(:));
    
    % 得到grad的直方图,并使用高的百分比(99.9%)估计梯度的阈值
    h =imhist(grad);
    Q = percentile2i(h,0.999);
    % 用Q对梯度做阈值处理,形成标记图像,并且从f中提取梯度值比Q大的点,得到结果的直方图
    markerImage = grad>Q;
    subplot(233),imshow(markerImage),title('以99.9%进行阈值处理后的梯度幅值图像');
    fp = f.*markerImage;
    subplot(234),imshow(fp),title('原图像与梯度幅值乘积的图像');
    hp = imhist(fp);
    % 用结果的直方图得到Otsu阈值
    hp(1) = 0;
    subplot(235),bar(hp),title('原图像与梯度幅值乘积的直方图');
    T = otsuthresh(hp);
    T*(numel(hp)-1)
    g = imbinarize(f,T);
    subplot(236),imshow(g),title('改进边缘后的图像')
    

    在这里插入图片描述
    示例2:用拉普拉斯边缘信息改进全局阈值处理

    f = im2double(rgb2gray(imread('test_img.png')));
    subplot(231),imshow(f,[]),title('原始图像');
    subplot(232),imhist(f),title('原始图像的直方图');
    hf = imhist(f);
    Tf = graythresh(f);
    gf = imbinarize(f,Tf);
    subplot(233),imshow(gf,[]),title('对原始图像进行分割的图像');
    w = [-1 -1 -1;-1 8 -1;-1 -1 -1];
    lap = abs(imfilter(f,w,'replicate'));
    lap = lap/max(lap(:));
    h = imhist(lap);
    Q = percentile2i(h,0.995);
    markerImage = lap >Q;
    fp = f.*markerImage;
    subplot(234),imshow(fp,[]),title('标记图像与原图像的乘积');
    hp = imhist(fp);
    hp(1) =0;
    subplot(235),bar(hp)
    T = otsuthresh(hp);
    g = imbinarize(f,T);
    subplot(236),imshow(g,[]),title('修改边缘后的阈值处理图像')
    

    在这里插入图片描述

    3.6 基于局部统计的可变阈值处理

    背景光照不均匀的情况下,或者在有多个主要物体灰度的情况下,此时使用全局阈值处理则不太合适,为此引入局部可变阈值处理
    我们用一幅图像中每个点(x,y)的邻域中像素的标准差和均值,分别表示局部对比度和平均灰度这两个描述子。
    令σxy和mxy表示一幅图像中以坐标(x,y)为中心的邻域Sxy所包含的像素集合的标准差和均值。则可变局部阈值Txy为:
    在这里插入图片描述
    其中a和b为非负常数,a和b的值通常是通过实验确定的。且
    在这里插入图片描述
    其中mG为全局图像均值。

    实现
    计算局部标准差,使用函数stdfilt进行操作

    g = stdfilt(f,nhood)
    nhood 是由 0 和 1 组成的数组,其中非零元素指定用于计算局部标准差所用的领域。nhood 的尺寸在每个维度上必须是奇数,默认值是 ones(3)。
    

    计算局部均值,可以用函数localmean进行操作
    localmean.m

    function mean = localmean(f,nhood)
    if nargin ==1
        nhood = ones(3)/9;
    else
        nhood = nhood/sum(nhood(:));
    end
    mean = imfilter(tofloat(f),nhood,'replicate');
    

    局部阈值处理操作,用函数localthresh来操作
    localthresh.m

    function g = localthresh(f,nhood,a,b,meantype)
    f = tofloat(f);
    SIG = stdfilt(f,nhood);
    if nargin == 5 && strcmp(meantype,'global')
        MEAN = mean2(f);
    else
        MEAN = localmean(f,nhood);
    end
    g = (f > a*SIG) & (f > b*MEAN);
    

    示例

    f = im2double(rgb2gray(imread('test_img.png')));
    subplot(2,2,1),imshow(f);title('(a) 原图');
    [TGlobal] = graythresh(f); 
    gGlobal = imbinarize(f, TGlobal); 
    subplot(2,2,2),imshow(gGlobal);title('(b)用 Otsus 方法分割的图像');
    g = localthresh(f, ones(3), 30, 1.5, 'global'); 
    SIG = stdfilt(f, ones(3));
    subplot(2,2,3), imshow(SIG, [ ]) ;title('(c) 局部标准差图像');
    subplot(2,2,4),imshow(g);title('(d)  用局部阈值处理分割的图像 ');
    

    在这里插入图片描述
    当背景接近于常数,并且所有物体的灰度高于或低于背景灰度时,选择全局均值一般会得到较好的结果。

    3.7 使用移动平均的图像阈值处理

    移动平均是一种局部阈值处理方法,该方法以一幅图像的扫描行计算移动平均值为基础。移动平均分割能减少光照偏差,且计算简单。当感兴趣的物体与图像尺寸相比较小(或较细)时,基于移动平均的阈值处理会工作的很好。适合于处理打印图像和手写文本图像。
    为减少光照偏差,扫描是以Z字形模式逐行执行的。令zk+1表示在扫描序列中第k+1步遇到的点的灰度。该新点处的移动平均(平均灰度)如下:
    在这里插入图片描述
    n表示用于计算平均的点的数量。
    使用下式进行分割:
    在这里插入图片描述
    其中K是[0,1]区间的常数, mxy是输入图像中点(x,y)处的移动平均(即上式的计算平均移动)。
    使用movingthresh函数来实现移动平均操作:
    movingthresh.m

    function g = movingthresh(f, n, K)
    f = im2double(f);
    [M, N] = size(f);
    if (n < 1) || (rem(n, 1) ~= 0)
    error('n must be an integer >= 1.')
    end
    if K < 0 || K > 1
    error('K must be a fraction in the range [0, 1].')
    end
    f(2:2:end, :) = fliplr(f(2:2:end, :));
    f = f'; 	% Still a matrix.
    f = f(:)'; % Convert to row vector for use in function filter.
    maf = ones(1, n)/n; 	% The 1-D moving average filter.
    ma = filter(maf, 1, f); % Computation of moving average.
    g = f > K * ma;
    g = reshape(g, N, M)';
    g(2:2:end, :) = fliplr(g(2:2:end, :));
    

    示例:用移动平均的文档阈值处理

    f = imread('noisy-text.png');
    f = rgb2gray(f);
    % 使用Ostu全局阈值分割
    T =graythresh(f);
    g1 = imbinarize(f,T);
    % 移动平均局部阈值分割
    g2 = movingthresh(f,20,0.5);
    figure;
    subplot(1,3,1); imshow(f); title('原始图像');
    subplot(1,3,2); imshow(g1); title('用otsu全局阈值分割后的图像');
    subplot(1,3,3); imshow(g2); title('用移动平均局部阈值分割后的图像');
    

    在这里插入图片描述
    关于n的选取,经验之一是令n等于平均笔画宽度的5倍,在此图像中,平均宽度为4像素,故n取20。
    原图是一幅以斑点灰度模式遮蔽的手写文本图像,用Ostu方法处理不能克服灰度的变化。很明显使用移动平均局部阈值处理分割能得到更好的效果。

    四、基于区域的分割

    以直接寻找区域为基础的分割技术。

    4.1 区域生长

    区域生长:根据预先定义的生长准则,将像素或子区域组合为更大区域的过程。
    基本思想:从一组“种子”点(种子点可以是单个像素,也可以为某个小区域)开始,将与种子性质相似的邻域像素或区域与种子点合并,形成新的种子点,重复此过程直到不能生长为止。种子点和邻域区域的相似性判据可以是灰度值、纹理、颜色等多种图像信息。
    算法步骤:
    1、选择合适的种子点;
    2、确定相似性准则即生长准则;
    3、确定生长停止条件。
    一般来说,在无像素或者区域满足加入生长区域的条件时,区域生长就会停止。

    用一个简单的例子来说明区域生长的过程:图(a)为原图像,数字表示像素的灰度值,以灰度为8的像素为初始的种子点。在8邻域内,生长准则是待测点灰度值和种子点灰度值相差为1或0。当合并到图d后,其他像素点已经不满足生长生长准则,故停止生长。
    在这里插入图片描述
    实际中,区域生长时经常还要考虑到生长的“历史”,还要根据区域的尺寸、形状等图像的全局性质来决定区域的合并。

    实现:
    通过自定义regiongrow函数来实现区域增长:

    [g,NR,SI,TI] = regiongrow(f,S,T)
    

    f:将要被分割的图像
    S: 可以是一个数组(与 f 大小相同)或一个标量。若 S 是一个数组,则它在所有种子点的坐标处必须是 1,而在其他位置为 0。若S 是一个标量,则它定义一个灰度值, f 中具有该值的所有点都会成为种子点。
    T: 可以是一个数(与f大小相同)组或一个标量。若 T 是数组,则它对 f 中的每个位置包含一个阈值。若 T 是标量,则它定义一个全局阈值。阈值用来测试图像中的像素是否与该种子或8连接种子足够相似
    例如:若S = a,且T = b,则我们将比较亮度,若 |像素亮度 - a| ≤ b,则称该像素类似于a。另外,若问题中的像素是8连接到一个或多个种子值的,则这个像素可看作是一个或多个区域的成员。
    g: 是分割后的图像,每个区域的成员都用整数标出。
    NR: 是所找到的区域的个数。
    SI: 是包含种子点的一幅图像。
    TI: 是一幅图像,该图像中包含在经过连通性处理前通过阈值测试的像素,即具有灰度zi且满足|zi - S| ≤T的点。
    SI和TI的大小均与f相同

    regiongrow.m

    function [g,NR,SI,TI]=regiongrow(f,S,T)
    %regiongrow执行区域生长
    f=double(f);
    %如果S是标量,则包含种子图像。
    if numel(S)==1
        SI=f==S;
        S1=S;
    else
        %S是一个数组。排除重复,在以下编码部分,连接种子位置去减少循环执行数量。
        SI=bwmorph(3,'shrink',Inf);
        J=find(SI);
        S1 = f(J);   % 种子数组的值
    end
    
    TI=false(size(f));
    for K=1:length(S1)
        seedvalue=S1(K);
        S=abs(f-seedvalue)<=T;
        TI=TI|S;
    end
    %使用SI的函数重构作为标记图像去获得与S中每个种子对应的区域
    %函数bwlabel给每个连接区域分配不同的整数
    [g,NR]=bwlabel(imreconstruct(SI,TI));
    

    示例:区域生长对焊接孔隙检测
    初始种子点:已知焊缝缺陷区域中的一些像素往往有最大的数字值(该情况下是255)

    f = imread('defective_weld.tif');
    [g,NR,SI,TI]=regiongrow(f,255,65);  %种子的像素值为255,65为阈值
    figure;
    subplot(2,2,1); imshow(f); title('原始图像');
    subplot(2,2,2); imshow(SI); title('种子点的图像');
    subplot(2,2,3); imshow(TI); title('所有通过阈值测试的像素');
    subplot(2,2,4); imshow(g); title('对所有连接到种子点的像素进行8连通分析后的结果');
    

    在这里插入图片描述

    4.2 区域分裂与聚合

    将一幅图像细分为一组任意的不相交区域,然后聚合和/或分裂这些区域。
    令R表示整幅图像区域,选择一个属性Q,对R进行分割就是依次将它细分为越来越小的四象限区域,以便任何区域Ri都有Q(Ri) = TRUE。
    从整个区域开始,若Q(R) = FALSE,则将其分割为4个象限区域,若分割后象限区域Q依旧为FALSE,则将对应的象限再次细分为四个子象限区域,以此类推。分割的过程可以用一个四叉树形式直观的表示。
    细分完成后,对满足属性Q的组合像素的邻接区域进行聚合,即Q(Rj∪Rk) = TRUE时,对两区域进行聚合。
    在这里插入图片描述
    算法小结:
    1、把满足Q(Ri) = FALSE的任何区域Ri分裂为4个不相交的象限区域;
    2、不可能进一步分裂时,对满足条件Q(Rj∪Rk) = TRUE的任意两个邻接区域Rj和Rk进行聚;
    3、无法进一步聚合时,停止操作。
    习惯上要规定一个不能再进一步执行分裂的最小四象限区的尺寸。

    实现
    使用函数qtdecomp 进行四叉树分解的操作

    S = qtdecomp(I,threshold,[mindim maxdim])
    I:输入的灰度图像
    threshold:分割成的子块中允许的阈值,默认为0。如果子块中最大元素和最小元素的差值小于该阈值,就认为该子块在属性Q上为TRUE,即不再进行进一步细分。对double型矩阵,threshold将直接作为阈值,而对uint8和uint16类型的矩阵,threshold将被乘以255和65535以作为实际的阈值。对图像而言,threshold的取值范围时0-1。
    [mindim maxdim]:尺度阈值。mindim可屏蔽函数对尺度上小于mindim的子块的处理,不管该子块在属性Q上是否为真;若参数为[mindim maxdim],则表示不产生小于mindim尺度的子块,也不保留大于mixdim尺度的子块。maxdim/mindim必须是2的整数次幂
    S:一个稀疏矩阵,在每个子块的左上角给出了子块的大小
    

    用一个小例子来说明一下四叉树分解结果矩阵S
    在这里插入图片描述
    红线分割开的就是对原图像分裂得到的子块,绿色框出每个块中左上角的非零元素就是块的大小。
    注: qtdecomp函数主要适用于边长是2的整数次幂的正方形图像;对长宽不是2的整数次幂的图像,分解可能无法进行到底。
    为了在四叉树分解后得到指定大小子块的像素及位置信息,使用函数qtgetblk

    [vals,r,c] = qtgetblk(I,S,dim)  
    I:输入的灰度图像
    z:由qtdecomp返回的稀疏矩阵
    dim:指定的子块大小
    vals:dim*dim*k的三维矩阵,包含I中所有符合条件的子块数据。其中k为符合条件的dim*dim大小的子块的个数,vals(:,:,i)表示符合条件的第i个子块的具体内容。
    r和c:均为列向量,分别表示图像I中符合条件子块左上角的行索引和列索引
    

    同样,用一个小例子来直观的感受一下

    >> I = [1    1    1    1    2    3    6    6
         1    1    2    1    4    5    6    8
         1    1    1    1   10   15    7    7
         1    1    1    1   20   25    7    7
        20   22   20   22    1    2    3    4
        20   22   22   20    5    6    7    8
        20   22   20   20    9   10   11   12
        22   22   20   20   13   14   15   16];
    >> S = qtdecomp(I,5);  % 对double类型的数据,5可以直接作为阈值
    >> [vals,r,c] = qtgetblk(I,S,4);  %  得到所有符合条件的大小为4*4的子块数据
    >> k = size(vals,3)   %  k为符合条件的子块的个数
    k =
         2
    >> % 显示第一个子块的数据
    >> vals(:,:,1)
    ans =
         1     1     1     1
         1     1     2     1
         1     1     1     1
         1     1     1     1
    >> % 显示第二个子块的数据
    >> vals(:,:,2)
    ans =
        20    22    20    22
        20    22    22    20
        20    22    20    20
        22    22    20    20
    

    使用qtsetblk函数将四叉树分解所得到的的子块中符合条件的部分全部替换为指定的子块
    参考:qtsetblk函数

    J = qtsetblk(I,S,dim,vals)
    I:输入的灰度图像
    S:I经过qtdecomp函数返回的稀疏矩阵
    dim:指定的子块大小
    vals:dim*dim*k的三维矩阵,包含了用来替换原有子块的新子块信息。k为图像I中大小为dim*dim的子块的总数,vals(:,:,i)表示要替换的第i个子块。
    J:经过子块替换的新图像
    

    用小例子进行直观感受:

    >> I = [1    1    1    1    2    3    6    6
         1    1    2    1    4    5    6    8
         1    1    1    1   10   15    7    7
         1    1    1    1   20   25    7    7
        20   22   20   22    1    2    3    4
        20   22   22   20    5    6    7    8
        20   22   20   20    9   10   11   12
        22   22   20   20   13   14   15   16];
    >> S = qtdecomp(I,5);
    >> vals = qtgetblk(I,S,4);
    >> newvals = cat(3,zeros(4),ones(4));  % cai将矩阵合成三维,第一层为zeros(4),第二层为ones(4)
    >> J = qtsetblk(I,S,4,newvals)
    J =
         0     0     0     0     2     3     6     6
         0     0     0     0     4     5     6     8
         0     0     0     0    10    15     7     7
         0     0     0     0    20    25     7     7
         1     1     1     1     1     2     3     4
         1     1     1     1     5     6     7     8
         1     1     1     1     9    10    11    12
         1     1     1     1    13    14    15    16
    

    由结果可以明显的看出,对原图像中大小为4*4的所有子块进行的替换。
    注: size(newvals,3)必须和I中指定大小的子块的数目相同,否则会出错。

    示例:对图像进行四叉树分解,并以图像形式显示得到的稀疏矩阵。

    I = imread('rice.tif');
    I1 = imresize(I,[512,512]);  % 改变原图的大小
    S = qtdecomp(I1,0.1);
    S1 = full(S);  % 将稀疏矩阵转换为普通矩阵
    
    figure;
    subplot(1,3,1); imshow(I); title('原图')
    subplot(1,3,2); imshow(I1); title('改变原图的大小');
    subplot(1,3,3); imshow(S1); title('四叉树分解结果');
    

    由于原图的大小为257257,直接对原图使用qtdecomp函数进行分解会报错,如下:
    在这里插入图片描述
    尝试将原图的大小修改为2的幂次,如代码中的512
    512,经测试,将阈值设置为0.1能得到较好的分割效果。
    在这里插入图片描述

    五、用形态学分水岭的分割

    分水岭算法(watershed)是一种借鉴了形态学理论的分割方法,在该方法中,将一幅图像看成一个拓扑地形图,其中图像的灰度值对应地形高度值。高灰度值对应山峰,低灰度值对应山谷。水总是朝地势低的地方流动,直到某一局部低洼处才停下来,这个低洼处被称为汇水盆地。最终所有的水会分聚在不同的汇水盆地,汇水盆地之间的山脊被称为分水岭。对图像进行分割,就是要在灰度图像中找出不同的“汇水盆地”和“分水岭”。也就是感兴趣的区域内部及其边缘。
    在这里插入图片描述
    一般考虑到各区域内部像素的灰度比较接近,而相邻区域像素间的灰度差距较大,可以先计算一幅图像的梯度图,再寻找梯度图的分水岭。在梯度图中,小梯度值对应区域内部,大梯度值对应区域的边界,分水岭算法就是在寻找大梯度值像素的位置。
    在这里插入图片描述
    直接应用分水岭分割算法的效果往往并不好,通常会由于噪声和梯度的其他局部不规则性造成过度分割。更有甚者,过度分割可能导致不可用的结果。一种解决方案是,通过融入预处理步骤来限制允许存在的区域数量。
    用于控制过度分割的一种方法是基于标记的。标记是指属于一幅图像的连通分量。与感兴趣物体相联系的标记称为内部标记,与背景相关联的标记称为外部标记。

    实现

    有关计算梯度图像的简单说明:以sobel算子为例

    f = imread('gel-image.tif');
    h=fspecial('sobel');
    fd=double(f);
    g=sqrt(imfilter(fd,h,'replicate').^2+imfilter(fd,h','replicate').^2);
    g1 = imgradient(f); % 使用imgradient函数也可直接得到图像的梯度
    g2 = imcomplement(imsubtract(g,g1));
    figure;
    subplot(2,2,1); imshow(f); title('原图');
    subplot(2,2,2); imshow(g,[]); title('用sobel算子计算图像梯度')
    subplot(2,2,3); imshow(g1,[]); title('用imgradient函数计算图像梯度')
    subplot(2,2,4); imshow(g2,[]); title('两种计算梯度方法的差值');
    

    在这里插入图片描述
    1、计算梯度可以使用sobel算子,然后分别求水平方向和竖直方向的偏导数,再计算梯度;也可以使用imgradient函数直接得到图像的梯度。
    2、用sobel算子计算梯度时,首先要将原图像转换为double类型,因为sqrt函数不支持uint8类型(用imread函数读入的图像为uint8类型);imgradient函数可直接对uint8类型图像进行处理,也可处理double类型图像
    3、将用两种方法得到的梯度图像相减,从图中可以看出效果是一样的。

    参考:标记控制的分水岭分割
    watershed函数——对图像进行分水岭算法处理,调用格式如下:

    L = watershed(A)
    或  L = watershed(I,conn)
    I:输入的需要处理的原图像
    conn:可选参数,用于指定分水岭算法将要考虑的邻域数量。对于二维图像,这个值可以为4或者8,该函数也可以处理三维图像。
    L:输出的矩阵中数值为0的元素表明不属于任何一个划分出的区域,而各个分水岭分割出的区域被用不同的序号表示,输出的矩阵时double类型 。
    

    imregionalmin函数——计算图像中大量局部最小区域的位置,调用格式如下:

    rm = imregionalmin(f)
    f:灰度图像
    rm:二值图像,二值图像的前景像素标记了局部最小区域的位置。
    

    imextendedmin函数——扩展最小值变换,即计算图像中比周围点更深的点的集合(通过某个高度阈值)。调用格式如下:

    im = imextendedmin(f,h)
    f:灰度图像
    h:高度阈值
    im:是一个二值图像,该二值图像的前景像素编辑了深局部最小值区域的位置。
    

    示例:控制标记符的分水岭分割

    f = imread('gel-image.tif');
    h=fspecial('sobel');
    fd=double(f);
    g=sqrt(imfilter(fd,h,'replicate').^2+imfilter(fd,h','replicate').^2);
    L=watershed(g);
    wr=L==0;  % L的0值像素是分水岭脊像素
    
    % 基于梯度图像进行分水岭变换,会看到由大量局部最小区域导致的过度分割
    rm=imregionalmin(g); % 计算图像中所有的局部小区域的位置
    im=imextendedmin(f,2); % 扩展极小值变换,用于计算比周围点更暗的图像中“低点”的集合
    fim=f;
    fim(im)=175;  % 上两行将在原图像上以灰色气泡的形式叠加扩展的局部最小区域位置
    
    Lim=watershed(bwdist(im));  % 通过计算内部标记符图像im的距离变换分水岭可以得到外部标记
    em=Lim==0;
    g2=imimposemin(g,im | em); % 使用 imimposemin 来修正梯度幅值图像,使其局部最小值只出现在前景和背景标记像素上。
    L2=watershed(g2);
    f2=f;
    f2(L2==0)=255;
    
    figure;
    subplot(2,4,1); imshow(f); title('原图');
    subplot(2,4,2); imshow(wr,[]);title('(a)梯度幅度图像进行分水岭变换图像');
    subplot(2,4,3); imshow(rm,[]);title('(b)对梯度幅度图像的局部最小区域');
    subplot(2,4,4); imshow(fim,[]);title('(c)内部标记符');
    subplot(2,4,5); imshow(em,[]);title('(d)外部标记符');
    subplot(2,4,6); imshow(g2,[]);title('(e)修改后的梯度幅度值图像');
    subplot(2,4,7); imshow(f2),title('(f)最后分割的结果');
    

    在这里插入图片描述

    展开全文
  • 摘要: 作为信息的最重要载体,电子文档处理的研究引起人们极大的兴趣。在任何文档处理系统中,预处理极为重要,其效果好坏会严重地影响其它模块的...此外,图像二值化本身也是数字图像处理中重要的基本问题。 本...
     

    摘要:

    作为信息的最重要载体,电子文档处理的研究引起人们极大的兴趣。在任何文档处理系统中,预处理极为重要,其效果好坏会严重地影响其它模块的工作。特别是灰度图像二值化效果的好坏,对识别效果以及其后的一切处理都有相当大的影响。原因之一是,任何物理传感输入都是灰度图像,文档处理系统的大多数模块却仅仅处理二值图像,图像二值化是必不可少的。此外,图像二值化本身也是数字图像处理中重要的基本问题。

    本文首先介绍完整的文档处理系统以及其预处理模块,然后将重点放在二值化问题上,给出图像阈值化方法的综述,并对全局阈值化、局部动态阈值化等方法的优点和缺点给出评价。本文尤其针对灰度变化比较复杂的文档图像,提出了一种改进的动态阈值化算法,并以仿真实验与以往的全局方法进行了比较,证实所发展新方法的优势。本文对该新算法在具体应用中会遇到的问题也做了估计,并提出解决相应问题的基本原则。

     

    关键词:预处理,灰度图像,二值化,阈值


    Direct Local-thresholding Method for Gray-level Document Images

    Abstract:

    People have great interest in the study of document processing, as electronic document is the most important information carrier nowadays. In any document processing system, pre-processing plays an important role in that it affects all the later modules in the system. Among various pre-processing, gray-level image to binary- image conversion, or image binarization is a key. First, an image obtained with a sensor is essentially a grey-level image, but most modules in the processing system only deal with binary image, and binarization becomes a necessity. Second, binarization has it’s own basic research interest.

    In this thesis, we first introduce the full document processing system and its preprocessing part; then we give special attention to the problem of binarization, and make a survey of various existing image thresholding methodologies, including global methods and local dynamic methods, with evaluations of their excellences and shortcomings. Then, aiming at the case of document image with complex gray-level variation, we bring forward an improved dynamic thresholding algorithm, named Direct Local-thresholding Method, which is proved to be better than former global method by simulating experiments. Finally, we discuss the problems may encountered in application systems, and put forward basic principles for possible solutions.

    Key words: preprocess, gray-level image, binarization, threshold

     


    第一章 文档图像预处理概述

    1.1 文档图像处理系统介绍

    物质、能量、信息,是客观世界中的三大要素。而文档作为信息的载体,在社会生活中占有着十分重要的地位。通常,我们可以把存储在计算机中或纸上的一切具有阅读意义的信息甚至承载着信息的纸统称为文档。文档可以分为两大类:文字/符号和图形/图像。

    为了能够更加简便有效地使用和管理信息,自二十世纪六十年代以来,人们进行了大量关于文档处理系统的研究。文档可以通过文档处理系统进入计算机,从而使得人们能够方便地对它们进行存储、管理、传输。文档处理涉及到很多领域,主要有光学字符识别(OCR),文档分析和理解,文档数据库,自然语言理解等等。

    文档进入计算机,是通过文档处理系统将输入的文档图像识别为具有意义的文字和图形来实现的。文档处理的一般流程如下[1]

    1.1 文档处理流程

    可以看出,在一个文档处理系统中,预处理作为后续工作的基础,是一个相当重要的部分。图像预处理的好坏对识别效果有相当大的影响。下面我们将对预处理部分作一下简要介绍。

    1.2 文档图像预处理

    在一个实际的文档处理系统中,预处理过程包含许多方面[2],例如:将输入的彩色图像转换为灰度图像;将灰度图像做二值化处理;背景噪音的消除;图像倾斜检测与校正;版面拆解等等。这里,我们将对其中比较重要的几个方面进行总体上的阐述。

    1.2.1彩色图像转换为灰度图像,以及灰度图像的二值化

    一篇存储在纸上的文档,要输入计算机中,首先要通过扫描仪扫描成为数字图像。它可以是彩色图像,也可以是灰度图像,或者是二值图像,取决于扫描的具体过程。但是一般来说,文档处理系统要处理的是灰度图像,而且很多成熟的图像处理算法和工具包也是针对灰度图像的。并且在进入识别阶段时,识别引擎一般是针对二值图像的。因此,我们必须对输入的图像进行处理,将彩色图像转化为灰度图像,再对其进行二值化。

    图像的这种转换必然会引起信息丢失,因此采用什么样的算法能够最大限度地保留识别时必需的信息(如字符的连通性),去掉不必要的背景信息和噪声,并且执行时间在实际可接受的范围内,是人们一直以来努力研究的问题。

    彩色图像转换为灰度图像的原理如下:彩色图像使用一个三维矢量 (R,G, B) 来表示一个像素点,而灰度图像用一个灰度级(gray level)来表示。因此这种转化可以看作是从一个三维矢量到一个一维矢量的投影操作。通常,可以用一个线性变换来完成这一过程,如下式:

    这里,R(x,y),G(x,y), B(x,y) 分别为像素(x,y)点的R,G,B三个分量的值,l,k,m为预先确定的参数。P(x,y)为求得的灰度值。

    图像的二值化,实际上是图像阈值化问题的一种特殊情况。图像的阈值分割是指,按照灰度级将图像空间划分成为与现实景物相对应的一些有意义的区域。在各个区域内部,灰度级是均匀的,而相邻区域的灰度级是不同的,其间存在着边界。阈值分割比较困难,但在机器视觉、文字识别、生物医学图像分析、指纹与印草鉴定、光学条纹判读以及军事目标识别等领域应用极为广泛。一般地,我们可以选择多个阈值,把图像的整个灰度级范围划分为几段,称之为多阈值分割。如果我们只选择单个阈值,将图像的灰度级范围划分为前景(或称目标)与背景两部分,就称为图像的二值化。

    在我们所研究的文档图像处理这一特定领域中,图像的二值化一般是指将灰度图像转化为只包含黑、白两个灰度的二值图像。文档图像可以看作是由文字、背景、图片三类对象组成,而我们最关心的是文字区域二值化后的结果。

    二值化问题是本文主要讨论的问题。关于它在后续章节中有更详细的阐述。

    1.2.2图像方向的自动检测与倾斜校正

    文档处理系统要求所处理的图像是正的,或者倾斜角度已知,否则许多对图像的操作,例如投影分析,图像分割等等,就无法进行。显然,扫描进计算机的文档图像无法保证一定是正的,因此需要利用倾斜检测和校正的方法对其进行处理。

    经常采用的倾斜角度检测的方法[3]有:基于文字行的检测方法,投影轮廓分析方法,和Hough变换方法等等。

    基于文字行的检测可以用于已知文字行方向(水平或竖直)的文档图像。它利用了对图像中连通体的分析。我们知道,连通体是一个灰度值相同的像素的集合,这个集合中任意两个像素之间都是8-近邻关系。可以用包含连通体内像素的最小矩形来表示连通体,它描述了连通体的大小和位置信息。如果已知文字行方向,我们就可以将连通体合并成文字行,并用直线逼近。该直线的倾斜角即为文字行的倾斜角。对整幅图像的文字行作同样分析,选出出现频率最高的角度,即可作为图像的倾斜角。

    投影操作也是一种基本的图像处理方法。将图像按一定方向作投影,可以得到在该方向坐标轴上分布的波形,它描述了图像沿该方向上的黑像素分布情况。如果图像文字行是水平的,那么沿水平方向的投影波形将具有明显的波峰和波谷。基于这一点,我们可以在候选倾斜角度范围内转动图像,直至出现明显的波峰和波谷为止,这时得到的角度就是倾斜角度。

    Hough变换方法是在倾斜检测中最常使用的方法,它抗噪声干扰的能力强,并且不受图像间隙干扰。它的原理是将直线从图像空间映射到参数空间,如下式:

    这样就将原图像中的直线映射到参数空间的一点,而原图像中的一点则对应着参数空间的一条正弦曲线。图像空间中任意两点所对应的正弦曲线在参数空间将相交于一点:,进而, 通过的直线上的所有点所对应的正弦曲线在参数空间都将相交于这一点。基于这一点,我们将直线检测问题转化为参数空间寻找局部最大值的问题。以上是标准Hough变换的基本思想。这种方法适用于文字行方向预先未知的情况。   

    1.2.3版面结构理解

    文档图像几何结构的理解也称为版面拆解,它是文档图像分析中的一个重要问题。版面拆解的目的是生成一个描述文档图像的层次结构:几何结构。它将图像分割成为具有相同特性的区域,为下一步的区域识别和文字识别做好准备。

    版面拆解的方法一般有如下三种:

    1.      自底向上(Bottom-Up,也称为数据驱动)。这类方法利用图像的局部信息,通过逐步将具有相同属性的区域合并,得到对文档版面的拆解。该方法能处理不同版面的文档和具有一定倾斜的文档,但是一般比较慢。

    2.      自顶向下(Top-Down,也称模型驱动)。该方法从图像全局出发,逐步对图像进行分割,最后得到图像的几何结构。对Manhattan式的版面[1],该方法快速而且有效;但是对复杂文档效果欠佳。影响自顶向下方法有效性的因素包括文字行位置的随意性,区域形状的不规则性以及文档图像的倾斜等。

    3.      综合方法(Hybrid)则尽量综合上述两种方法的特点,使算法的性能和适应性都得到提高。

    由南开大学机器智能研究所开发的识别工具包RTK7.0中,版面结构理解模块采用了基于模型序列的方法[4]。这种方法首先通过大量的观察,建立一些基本的文档结构模型,如:多节结构,多栏结构,标题加多栏结构,非平衡栏结构,等等。另外还有一些由这些基本结构模型组合而成的常见结构模型。实践证明,这些结构模型虽然比较简单,但是能够很好地表示实际应用中出现的版面结构。在对三个彼此独立的样本图像集合的测试中,其覆盖系数达到92.90%

    这一章我们从总体上阐述了关于文档处理系统及其预处理部分。从下一章开始,我们将重点研究灰度图像的阈值化问题。


    参考文献

    [1]     刘秋元,“文档识别和统计分析”,南开大学博士论文,1997

    [2]     李昀,“OCR系统前处理的研究”,南开大学硕士论文,2000

    李庆中,“OCR中字符切割及文字区域抽取的研究”,南开大学硕士论文,2001

    潘武模,“模型序列方法与文档版面结构理解

    [5]      ”,南开大学博士论文,2001


    第二章 传统的图像阈值化方法介绍

    2.1 图像分割与图像阈值化概述

    图像分割(image segmentation)是图像处理中的主要问题,也是计算机视觉领域低层次视觉中的主要问题,同时也是一个经典难题。关于图像分割技术,从二十世纪七十年代起,就有许多研究人员和学者进行了大量研究,目前已经有相当多的研究成果。但是由于问题的重要性和困难性,到目前还不存在通用的方法,也不存在判断分割是否成功的客观标准。

    所谓图像分割[1],是指将图像中具有特殊涵义的不同区域区分开来,这些区域互不相交,每一个区域都满足特定区域的一致性。

    如果N为自然数集,(x , y)为数字图像的空间坐标,G= {0,1,, l-1}为代表灰度值的正整数集,则图像函数可以定义为如下的映射f (x , y)代表坐标为(x , y)的像素点的灰度值。

    将一幅图像进行分割[2],就是将图像划分为满足如下条件的子区域

    (a)    即所有子区域构成了整幅图像;

    (b)   是连通的区域;

    (c)    即任意两个子区域不存在公共元素;

    (d)   区域满足一定的均匀性条件。

    均匀性一般指同一区域内的像素点之间的灰度值差异较小,或灰度值的变化较缓慢。

    图像的阈值选取是图像分割的关键技术。使用阈值分割,就是用一个或几个阈值将图像的灰度分布直方图划分成几个类,认为灰度值在同一个灰度类内的像素属于同一个物体。通常,选择一个阈值将灰度范围划分成目标和背景两类,称为图像的二值化。

    tG作为阈值,B = { }为一对灰度级,。用灰度级 t 作为阈值对图像函数 f (x , y)进行二值化的结果是产生一个二值的图像函数

    总之,阈值选取方法就是基于某种标准来确定一个 t的值

    一般来说,现有的阈值选取技术可以分为全局的和局部的阈值选取方法。全局的阈值选取是指根据整幅图像确定一个阈值。局部阈值选取方法是指将图像划分为若干子图像,根据每个子图像确定相应的阈值。

    全局阈值选取方法对噪音比较敏感,因此应用中一般采用局部阈值选取方法。但是这二者并无本质的不同。

    2.2 全局阈值选取方法

    全局阈值选取方法又可以进一步分为点依赖的(point-dependent)和区域依赖的(region-dependent)两种[4]。如果 的确定仅仅依赖于每个点的灰度值,则这种阈值选取的方法就是点依赖的。如果 的确定依赖于邻域点的局部特性(如:局部灰度值的分布),则这种阈值选取的方法就是区域依赖的。点依赖的全局阈值方法有:直方图方法;最大类间方差法;最小误差法;矩量保持法;最大熵方法等等。区域依赖的全局阈值方法则有:直方图变换法;概率松弛方法;共生矩阵法;灰度分布统计方法等等。下面,我们将分别介绍其中比较重要的方法,尤其重点介绍最大熵方法。

    2.2.1直方图与直方图变换方法

    直方图方法[3] [4]是指直接从原图像的灰度分布直方图上确定阈值,包括p-分位数法、最频值法和直方图凹面分析法。

    p-分位数法(p-tile) [4]可以说是最早的一种阈值方法。这种方法假设图像包含黑的目标和亮的背景,目标所占区域的百分数已知。该方法设定的阈值使至少(100-p)%的像素点在二值图像中匹配为目标。例如,已知一页文件上印刷文字约占有全页纸的25%,那么所选阈值应使得灰度级小于阈值的像素数目约占总数的25%。显然,这种方法对于事先未知目标点数占有像素总数百分比的图像是不适宜的。

    最频值法(mode ) [4]也是一种很常用的简单方法。在较理想状态下,图像中的目标和背景非常清楚,灰度直方图呈现明显的双峰状。这时可以选取两峰之间的谷底(最小值)对应的灰度值为阈值。但是,实际图像的情况往往比较复杂,而且可能有很多噪音干扰,直方图参差不齐,很难确定极大值和极小值。另外,这种方法也不适用于直方图中双峰值差别很大,或双峰间的谷宽广而平坦的情况,以及单峰直方图的情况。因此这种方法常常与其他方法结合使用。

    直方图凹面分析[4]RosenfeldTorre1983年提出的,它适用于直方图中不存在谷的图像。它的原理是,由于直方图中“谷”和“肩”都对应凹面,可以在“肩”处寻找一个好的阈值。设HS表示在一组灰度值0,1,,l-1上定义的直方图,h(0)h(1),…,h(l-1)分别表示直方图中这些灰度值对应的高度,这里 ,有h(i) 0。这样HS可以看成一个二维区域。为了确定HS的凹面,构造包含HS的最小凸多边形,由集合差HS确定HS的凹面。若h(i)(i)分别表示HS在灰度值i之处的高度,则当(i)h(i)取局部极大值时所对应的灰度值i可以考虑作为阈值。为了去除伪凹面(由噪声脉冲引入)引起的极大值,引入一个平衡度量

    .

    度量了直方图在灰度值i处的平衡性。对于伪凹面,它通常位于直方图的某一侧,其值比较小。因此,略去值较小时(i)h(i)的极大值,其余极大值所对应的灰度值可以选作阈值。这不一定是最优的,可以考虑在这些极大值附近的其他灰度值。

    上述几种方法都是直接从原始灰度分布直方图上获取阈值的。实际上,图像的原始灰度直方图往往十分粗糙,存在很多毛刺,如果直接使用直方图方法,常常不能得到理想的结果。基于此原因,人们提出了直方图变换法,20世纪70年代关于阈值选取的大部分工作都集中在这一方面。

    直方图变换法[3]并不直接选取阈值,而是将原始灰度直方图经过某种变换,使之具有较尖锐的峰和较深的谷,这样就可以应用前面所述的直方图方法选取阈值。直方图变换法的共同特点是根据像素点的局部特性,为图像中的像素加权,从而获得新的直方图。它们都假设图像可以看作由目标和背景两部分构成,每部分都具有一致的灰度分布。主要的直方图变换法有以下几种:[4]

    1)        PandaRosenfeld提出的仅由边缘值(灰度变化率)较低的像素点构成灰度直方图的方法。由于目标和背景内部像素点的边缘值一般较低,而目标与背景边界上像素点的边缘值较高,因此,由此构成的直方图与原始直方图相比,双峰基本保持不变,而谷将变得较深。Katz也曾指出过由边缘值较高的像素点构成的灰度直方图应呈现出一个单峰,该单峰对应的灰度级可选作阈值。

    2)        Mason等人提出的使用边缘检测算子(例如Laplacian算子、Robert算子)对直方图进行加权。均匀区域处像素的边缘值较低,给予较大的权;而边界领域中边缘值较低,给予较小的权。权因子以计算,这里为给定像素处的边缘值。

    3)        WezkaRosenfeld1979年提出,通过构造图像的灰度-边缘值二维散射图,并计算在灰度级轴上的各种加权投影,可以得到变换后的直方图。

    4)        Wu等人提出的四元树方法。它是基于灰度标准差较高的非均匀区域可划分为若干标准差较小的均匀区域这一事实。首先从原始图像开始,如果标准差超过一预定值,就将图像一分为四,然后对每个子图像重复这一过程,直至将原始图像划分成若干标准差较小的区域块。用每块的平均灰度级来近似代替该块的灰度值,由此得到的图像称之为Q-图像,其灰度直方图将峰谷分明。由于 Q-图像中较小的图像块势必位于区域边界的附近,且图像块较小,离边界越近。因此,进一步深化直方图中谷的一种方式是抑制这些较小的图像块,从而使阈值选择变得更为容易

    2.2.2最大类间方差法(ostu方法)

    最大类间方差方法[4]是由Ostu1978年提出的,它是一种基于判别式分析的方法。把图像中的像素按灰度级用阈值t划分成两类,即= {0,1,t}= {t+1,t+2,l-1}。若用分别表示类内,类间和总体方差,则通过使下列关于t的等价的判决准则函数达到最大来确定最佳阈值

    ,    ,    .

    三个准则函数中,最为简便,又因已知,与t值无关,因此最优阈值

    这里,

    ,    ,    .

    Ostu方法计算简单,稳定有效,是实际应用中经常采用的方法。但是它对噪声和目标大小十分敏感,仅对类间方差是单峰的图像有较好效果。当目标与背景的大小比例悬殊时,类间方差准则函数可能为双峰或多峰,此时Ostu方法就会失效。[5]

    2.2.3最小误差法

    最小误差法[3]中,灰度直方图视为目标与背景像素灰度构成的混合集概率密度函数p(g)的估计,通常假设混合集的两个分量p(g / i) (i = 1 , 2)服从均值为,标准差为的正态分布,先验概率为,即

    ,

    这里

    .

    求解下列二次方程可得到Bayes最小误差阈值

    .

    然而,参数,通常是未知的,为了克服估计这些参数的困难,KittlerIllingworth引入了如下的准则函数:

    这里

     ,      ,

     ,             ,

     ,      .

    最佳阈值,当J(t)取最小值时获得,即

    .

    最小误差法计算简便,受目标大小和噪声影响小,对于小目标图像,例如目标与背景大小之比低于1的图像也有很好的效果。该方法可用迭代方式实现,并可推广到多阈值的选取。

    2.2.4概率松弛法(Relaxation

    概率松弛法[6]假设图像是由“亮”的背景与“暗”的目标(或者相反)构成的,这样就可以使用一个阈值来分割。它对各个像素点及它们之间的空间相关性用概率加以描述,通过多次迭代使对象与背景区域得到明显区分,能够有效消除背景噪声。

    概率松弛法的基本思想是,首先把像素根据其灰度级按概率分成两类;然后,按照邻近像素的概率调整每个像素的概率,调整过程迭代进行,以使对于属于亮(或暗)区域的像素,(或暗)概率变得非常高。

    关于像素的初始分类,RosenfeldSmith提出了如下方法。设dl分别是最暗和最亮的灰度级,是像素的灰度级,则对,令

    这个方法虽然简单,但仅限于目标与背景灰度级分别位于灰度级直方图前后两半的情况。为此,他们又提出了另一个初始化策略。设m是平均灰度级,则

    Rosenfeld与S.W.Zucker 提出的概率松弛法迭代公式如下:

    设像素可以被初始化分成l类(这里l=2),则

    这里,是对第i个像素点的所有邻近点(h)计算下式:

     后求和取得的均值。是在第r次迭代中,对第i个像素点属于第j类的概率的估计。C(i,j;h,k)是相关系数,用于衡量事件(像素点ij类内,像素点hk类内)的相关性。其定义为:

    2.2.5最大熵方法

    信息论之父Shannon 在其经典著论《通信中的数学原理》中提出,一个系统的信息量,即熵,是对于系统的不确定性的度量。近年来,有很多人将信息理论应用于图像处理和模式识别领域。Pun[7]使用Shannon关于熵的概念,在假设图像仅完全由其灰度分布直方图表示的基础上,定义了图像的熵,并用这一量度来实现目标/背景的分离。Kapur等人[8]也使用熵来进行图像分割,不同的是,他们使用的不是整个直方图的单一概率分布,而是使用目标与背景两个相互独立的概率分布。Abutaleb[9]则考虑灰度的空间分布,即图像的更高阶熵,他使用像素灰度值及其邻近像素的平均灰度值来计算图像的二维直方图,并由二维概率分布得到二维熵。Pal等人根据Deluca Termini定义的模糊集上的非概率性熵度量,实现了图像阈值选取的算法。另外,Pal Bhandari[10]使用基于poisson分布的、前景与背景的条件熵,并证明在阈值选取上优于前述算法。

    最大熵方法对不同目标大小和信噪比(SNR)的图像均能产生很好的分割效果。它受目标大小的影响小,可以用于小目标分割。但是它运算量大,运算速度比较慢。[5]

    本文主要研究使用最大熵方法进行图像的阈值化。下面将会对某些经典算法进行介绍。

    1.1.1.1    2.2.5.1 Shannon关于熵的定义

    Shannon认为,任何随机事件都含有某种不确定性因素,称为不确定性(uncertainty)。当这一事件发生后,人们在这一过程中获得了其不确定性,即信息。但是不同的随机现象中包含的不确定性是不同的,因而信息的多少也是不同的。所以信息的度量是信息论的最基本的问题。

    一般地,定义一个具有n个状态的系统的熵如下:

    为系统的第i个状态

    发生的概率。

    PalPal[11]指出了上述定义的一些缺点,重新定义了熵的形式如下:

    因此,我们可以将熵的形式写作:

    这里依赖于所使用的定义。

    从熵的定义可以看出,系统越不稳定,熵就越大,也就是说概率分布一致性强的系统具有较小的熵。

    1.1.1.2    2.2.5.2 Pun的最大熵方法

    Pun[7]提出的方法如下:

    t为图像分割的阈值,定义后验熵:

    ,

    ,

    这里,可以分别看作二值化以后,与黑,白像素有关的后验信息的度量。

    在知道了灰度直方图的先验熵的情况下,Pun 提出一种确定最优阈值的算法,它通过使后验熵

    的上界取最大值而得到。Pun 指出使后验熵取最大值相当于使如下估价函数

    取最大值。这里

    , ,.

    另外,Pun 还提出一种使用各向异性系数来选取阈值的方法,这里

    这里m是使得

    的最小整数。选取最优阈值,使得

    .

    但是,Kapur等人发现,这种方法总是得到,从而引入了不必要的偏移。

    1.1.1.3    2.2.5.3 Kapur,Sahoo和Wong的方法

    Kapur等人[8]提出的这种方法中,由原始图像的灰级直方图定义了目标和背景两个相互独立的概率分布如下:

    ,

    这里,t为阈值,。定义:

    ,

    .

    求得阈值使熵达到最大

    .

    这时可将作为最佳阈值。

    1.1.1.4    2.2.5.4 Abutaleb的二维熵算法

    Abutaleb[9]使用图像的像素灰度值及其邻近像素的平均灰度值来计算图像的二维熵。令表示一幅M×N的图像,一维阈值方法是用灰度T来将图像中所有像素点划分为两类,二值化后的图像可以用下式表示:

    g(x,y)分别为局部平均灰度值函数和二维阈值函数。如果图像由灰度值T和局部平均灰度值S划分,那么

    2.1 使用阈值矢量(T,S)分割的二维灰度直方图

    2.1是二维的(灰度值,局部均值)直方图。原点定义于左上方,灰度值自左向右增长,局部均值自上向下增长。由于灰度值数目为L,所以直方图中共有个元素。每个元素都代表相应的事件对(灰度值,局部均值)的发生。如果图像用二维阈值(T,S)来分割,那么上述直方图将被划分为4个象限。由于目标或背景内部的像素分布是一致的,所以它们的分布主要在从原点到(L -1,L-1)的对角线附近,也就是说,0象限和1象限分别包含了目标类与背景类的像素分布。而2象限和3象限则主要包含了边缘处的像素点和噪音像素点的分布。

    二维熵阈值选取的基本思想是选取阈值矢量(T,S),使得目标类和背景类的后验熵可以达到最大值。由于这两类的分布是互相独立的,我们可以定义两类中每个事件对(灰度值,局部均值)的概率如下:

    两个类的熵H0(T,S)H1(T,S)定义如下:

    选择阈值向量(T,S),使得

    H(T,S)= H0(T,S)+ H1(T,S) 达到最大,这时的(T,S)即为最佳二维阈值。

    另外,Brink提出另一种优化的算法,通过使目标类和背景类的熵中较小者达到最大值来求得阈值矢量。二维熵方法在效果上要优于一维熵算法,但是其运算时间复杂度为,非常耗费时间。基于此,ChenWen等人在Brink的算法的基础上提出了一种快速二维熵算法。[12]

    这种快速算法的思想如下:二维熵阈值方法的时间消耗的关键在于计算(T,S)对的熵,那么可以估计出一个可能的阈值向量的集合,只对集合中的元素应用计算熵的操作。算法分为两段式操作:第一阶段,求出阈值向量的候选集合。把原始图像的灰度范围除以,来进行量化,这样,量化图像的每个灰度级代表原始图像的一个大小的连续灰度级。量化阈值向量代表(T,S)对的一个区域,这就是候选阈值集合。第二阶段,在这个区域里来计算全局阈值向量。实验证明,最优阈值总是落在候选集合内。这个算法的时间复杂度为,大大提高了算法执行速度。

    1.1.1.5    2.2.5.5 Pal和Bhandari的条件熵方法

    PalBhandari[10]1993年提出条件熵方法,并证明比Pun的方法、Kapur等人的方法、PalPal的方法更优。

    条件熵是以2-序共存矩阵(second-order co-occurrence matrix)的形式定义的。对于如前所述的图像f(x,y),它的共存矩阵是一个L×L维的矩阵,表示了邻近像素间灰度跳变的信息。也就是说,矩阵中的元素表示从灰度i跳变到灰度j的次数。其定义如下:

    a表示图像f(x,y)中的像素点(i,j)b表示其8个邻近像素之一,即:

    定义共存矩阵中的元素,这里

    条件熵的定义如下:

    考虑两个独立实验ABAm种可能性Bn种可能性。已知发生的条件下,A的条件熵为:

    这里,是已知发生的条件下,发生的条件概率。

    由此,以B为条件,A的条件熵为:

    这里是事件对的联合概率分布。

    表示在邻近像素点灰度值j属于背景类的条件下,灰度值i属于目标类的条件概率, 。这样,对于给定阈值T,在已知背景类的条件下,目标类的条件熵定义为

    这里

    是事件对(i,j)发生的频率,也即前述共存矩阵中的第i行第j列的元素。

    同样可定义已知目标类条件下,背景类的条件熵

    这样,分割后图像的总体条件熵为

    选择阈值T使得H达到最大值,此时的T即为最佳阈值。

    关于这一点,可简单证明如下:[13]

    2.2 用阈值分割共存矩阵

    th为目标/背景分割的正确阈值。用th来分割共存矩阵,如图2.2所示,显然,目标类和背景类内部的像素点分别主要分布在区域13内,而边缘像素点和噪音点主要分别在区域24内。由于目标类和背景类内部的像素点灰度分布是一致的,所以在同一类别内,从某一灰度值跳变到另一灰度值的频率会很高,而从一个类别内的灰度值跳变到另一类别的灰度值的频率会比较低。也即区域13内的频率比较高,区域24内的频率比较低,而且它们都具有一致分布。现在,令我们假定的阈值为T,且有T<th,那么用T来分割共存矩阵,就会使区域2的频率变高(实际上这些跳变是发生在目标类内部的),而原来区域2内的跳变频率,即从目标类到背景类的跳变频率,仍然会比较低,这样区域2就会具有不平衡的概率分布,从而使其熵值降低。同理,区域1会保持其概率分布的一致性,而区域34受到的影响也会降低其熵值,因此总体上的熵值必然降低。同样分析T>th的情况,也会得出熵值降低的结果。因此,H关于T的最大值是目标/背景分离的最佳阈值。

    2.3 局部阈值方法与动态阈值选取

    对于目标和背景比较清楚的图像,全局阈值化方法可以取得较好结果。但是如果图像的背景不均匀,或是目标灰度变化率比较大,全局方法一般就不适用了。如下图所示:

    2.3 a

    2.3 b: 2.3 a使用全局方法阈值化的结果

    2.4 a

    2.4 b:2.4a使用全局方法阈值化后的结果

    这种情况下,尽管局部看来,目标与背景是可分的,但是无法得到一个适用于整幅图像的全局阈值。因此,人们提出了很多动态的局部阈值化算法,也称自适应阈值化算法。所谓动态是指,根据每个像素及其邻域像素的灰度值情况动态地计算分割所需的阈值。常用的动态阈值化算法有ChowKaneko的方法[14]YanowitzBruckstein的方法[15],等等。由于动态阈值化方法常常需要对图像中每个像素点都计算阈值,也就是说,对整幅图像求出一个阈值面(通常是曲面),因此计算量很大,运算速度一般比较慢。

    下面将介绍一些常见的动态阈值化算法。

    2.3.1 Chow和Kaneko的方法

    ChowKaneko[14]1972年提出的方法是通过灰度分布的局部信息来计算整幅图像的自适应阈值面(adaptive threshold surface)。这种方法的基本思想如下:首先,将整幅图像划分成为许多互不重叠的、相等的小单元,组成规则的网格,对每个小单元计算其局部灰度直方图。然后,判断哪些直方图具有双峰,并且利用它们获得对应小单元的局部阈值,对这些局部阈值进行插值,得到整幅图像的阈值面。这种算法比起全局方法有很大的改进,但是其缺陷也很明显:因为如果小单元的面积取得太小的话,用于产生直方图的像素数目也会太少,其统计意义也随之变得不明显。另外由于网格的划分与图像的内容并无关系,由小单元得到的局部阈值比较随意,而不是在比较有意义的层次上。还有,由于可能受到噪音等影响,小单元可能会整个落在目标或背景区域内,而直方图仍呈双峰分布,这时获得的阈值面就会发生严重错误。

    2.3.2 Yanowitz和Bruckstein的方法

    1989年,YanowitzBruckstein[15]提出了一种不仅仅依靠灰度分布,而是通过求灰度梯度的最大值来确定阈值面的方法。该算法的主要步骤如下:

    1)        对图像做平滑操作,即用每个像素的邻域灰度均值来替代实际灰度值。平滑后的图像更容易区分目标/背景。

    2)        从平滑后的图像产生灰度梯度图像。

    3)        对梯度图像应用阈值化和基于局部最大值的细化操作,获取具有局部最大梯度值的像素点。这些点的原始灰度值通常可以作为局部阈值的候选值。

    4)        对平滑后的图像在具有局部梯度最大值的点出取样,然后对这些点在整个图像范围内进行插值,构成阈值面。

    5)        用此阈值面进行图像分割。

    2.3.3 Sauvola和Pietikäinen的方法

    SauvolaPietikäinen[16]2000年提出了一种针对于文档图像的二值化方法。这种方法先将文档页面中的内容分为背景,图片和文字几类,然后对不同类的部分调用不同的二值化算法。

    首先,将文档图像划分为大小相等的矩形窗口,可以取大约10-20像素宽。对每个小窗口计算像素平均灰度值和跳变差(transient difference)两个特征。跳变差是局部对比度变化的量度,它的值在[0,1]区间上。设小窗口宽度为n,则在小窗口内计算跳变差的公式如下:

    f(i,j)为图像f(x,y)(i,j)像素点的灰度值。

    跳变差特征可以被分为以下4种:

    跳变差值

    特征

    一致

    几乎一致

    变化

    跳变

    “一致”的和“几乎一致”的跳变差特征对应着背景和图片内容,而“变化”的和“跳变”的则对应着文本类内容。对于每个小窗口,其二值化算法选择的规则如下:

    1)        如果灰度均值比较高,且灰度直方图的全局峰值与之落在同一个四分之一直方图内,且跳变差特征为“跳变”,那么使用基于模糊决策的阈值化方法(SDM)

    2)        如果灰度均值中等,且灰度直方图的全局峰值与之不在同一个四分之一直方图内,且跳变差特征为“一致”,那么使用基于直方图分析的阈值化方法(TBM)

    2.5 SauvolaPietikäinen的自适应文档图像二值化方法流程

    在这一章里,我们主要介绍了已有的一些图像阈值化方法和一些基本概念。从下一章开始,我们将针对某类较特殊的文档图像,提出一种新的二值化算法。


    参考文献

    [1]     罗希平等,“图像分割方法综述”,模式识别与人工智能,1999vol.12,no.3, pp.301-312

    [2]     N.R.Pal and S.K.Pal, “A Review on Image Segmentation Techiniques”, Pattern Recognition, ,1993, vol. 26, no. 9, pp. 1,277-1,294

    [3]     吴全,朱兆达,“图像处理中灰度级阈值选取方法30年(1962-1992)的进展(一)”,数据采集与处理, 1993 8(3-4)pp.193-201

    [4]     P.K.Sahoo, S.Soltani, A.K.C. Wong, and Y.C.Chen, “A Survey of Thresholding Techniques”, Computer Graphics,Vision and Image Processing 1988 (41),pp.233-260

    [5]     刘文萍,吴立德,”图像分割中阈值选取方法比较研究”,模式识别与人工智能,1997vol.10,no.3,pp.271-277

    Rosenfeld,R.C.Smith,”Thresholding using relaxation”, IEEE Trans. Pattern Anla.Mach.Intell., PAMI-1981 3(5), pp.598-606

    [7]     Pun,T., “A new method for gray level picture thresholding using the entropy of the histogram”. Signal Processing ,1980 2(3),pp.223-237

    [8]     J.N.Kapur, P.K. Sahoo, A.K.C. Wong, “A new method for gray-level picture thresholding using the entropy of the histogram”, Computer Graphics,Vision and Image Processing , 1985,29, pp. 273-285

    [9]     Abutaleb, A.S., “Automatic thresholding of graylevel pictures using two-dimensional entropy”. Computer Graphics,Vision and Image Processing , 1989(47), pp.22-32

    [10]Pal,N.R.,Bhandari,D., “Image thresholding:Some new techniques”. Signal Processing , 1993,33(2),pp.139-158

    [11]Pal,N.R.,Pal,S.K., “Object background segmentation using new definition of entropy”. Proc. IEEE, 1989,vol.part E,pp.284-295

    [12]W.Chen,C.Wen,C.Yang, “A fast two-dimensional entropic thresholding algorithm”, Pattern Recognition , 1994, 27(7), pp.885-893

    [13]Sambhunath Biswas,Nikhil R. Pal, “On hierarchical segmentation for image compression”,Pattern Recognition Letters, 2000 (21),pp.131-144

    [14]C. K. Chow and T. Kaneko, "Automatic detection of the left ventricle from cineangiograms", Computers and Biomedical Research, , 1972,vol. 5, pp. 388--410

    [15] S.D.Yanowitz,A.M.Bruckstein, “A new method for image segmentation”, Computer Graphics,Vision and Image Processing , 1989,46, pp. 82-95

    [16] J.Sauvola,M.Pietikäinen, “Adaptive document image binarization”, Pattern Recognition ,2000(33), pp. 225-236


    第三章 改进的算法:直接局域二值化方法

    3.1 问题的提出

    如前所述,二值化是文档图像处理系统中预处理模块的一个重要部分。文档图像作为数字图像的一类,有其自身的特点。一般来说,文档图像的二值化中,目标/背景的分离也就是文字/背景的分离。在进行文档图像处理时,我们最关心的是恢复其文本内容,因此在进行二值化预处理时,最重要的也是保证文字区域的二值化效果。

    实际应用中,由于文档图像处理系统的实时性要求很高,为了算法的执行速度考虑,在预处理时尽量采用简单的二值化算法,一般使用一个全局阈值化方法。实践证明,对于大部分的报刊、杂志等非彩色样张,其二值化效果可以满足文字识别模块的需要。但是随着排版和印刷技术的发展,文档样张的形式越来越丰富多彩,文档中的背景、文字都不再是单一颜色,而常常具有各种纹理和多变的色彩。这对于我们的文档图像处理系统是很大的难题和挑战,因为这种情况下,原来的文档图像二值化处理显然已经不能满足需要,必须加以改进。

    3.1 a

    3.1 b

    3.1 c

    3.1是一些文字和背景颜色不均匀的文档样张例子,这些样张在目前的杂志、学术期刊、论文、报纸上都很常见。一般地,颜色使用较复杂的样张可以大致分为如下几类:

    1)    正文部分文字颜色不一致。

    2)    背景颜色不一致,呈现出分块状,但在每一块内颜色一致。

    3)    背景有纹理或图案。

    可以看出,图3.1 a12两种情况的综合,图3.1 b属于第3种情况,图3.1 c则属于第2种情况。

    本文所要解决的问题主要针对于类似于图3.1 c的情况,即正文文本颜色一致,但是背景颜色有块状变化的文档图像的二值化问题。这种样张常见于各种学术期刊和论文中。

    3.2 对问题的分析及解决方案的提出

    考虑到文档图像处理系统的要求及本类样张的特点,本文在对各种图像阈值化算法进行综合比较、分析后,提出了一种改进的动态阈值化方法,我们把它叫做直接局域二值化方法。

    首先,基于文档处理系统对实时性的要求,我们不能采用运算量过大的复杂算法。如果对图像中每一像素点都进行动态的阈值计算,势必大大影响计算速度,无法应用于实际系统中。因此本文认为,对图像采取量化的阈值计算比较合适,即将图像划分为许多小区域,对每个小区域分别应用阈值化算法选取一个阈值,再从这些局部阈值得到整幅图像的阈值面。

    其次,对每个小区域的阈值化算法,可以采用最大熵方法。这是因为文档图像中目标(文本)所占比例较小,而且经常有一些噪音干扰,而最大熵方法对于目标大小和噪音不敏感,可以取得较好的分割效果。但是最大熵方法计算量大,计算速度慢,因此有必要做一些优化处理。本文提出的方法采用了PalBhandari[1]的条件熵方法(参见第二章2.2.5.4),使用灰度的空间分布及像素的空间相关性进行阈值分割,并对其进行了一些改进,使其运算时间在可以接受的范围内。

    下面将对本文提出的方法进行介绍。

    3.3 直接局域二值化方法

    本算法主要由以下步骤组成:

    1.      将原始图像划分为互不覆盖的、分别具有一定一致性的子区域;

    2.      对各个子区域动态计算阈值;

    3.      用得到的阈值点构成整幅图像的阈值面,并用此阈值面进行二值化操作。

    以下对各个步骤分别详细阐述。

    3.3.1将图像划分成为子区域

    将图像划分为小区域,是本算法的第一步,也是非常重要的一步。

    如何才能将图像大致划分成为意义一致的子区域,而计算代价又不致过大?回顾ChowKaneko[2]提出的最早的动态阈值化方法(参见本文第二章2.3.1),他们首先将图像划分为规则的、等大小的网格。这是最简单的方法,但是这种划分与图像的内容完全无关,导致划分后的网格意义不大。

    本文提出的区域划分方法是,先将图像用上述规则网格划分后,再根据各个网格计算出特征,对其分别进行进一步的合并、拆分等处理,得到基本意义一致的子区域。其主要步骤如下:

    1.          对于M×N的图像f(x,y),使用n×n大小的网格对其进行划分。这里n不能取得过大,不然就失去了划分的意义;而如果取得过小,计算特征时其统计意义不明显,而且会给后面的合并、拆分处理带来困难。根据观察,对于分辨率为300dpi,普通A4页面大小的常见文档图像,n100200像素宽时有较好效果。

    2.          对每个网格,计算下列特征:灰度均值,灰度方差,灰度跳变差[3]。公式如下:

    3.          对于每个网格,根据上述特征判断网格性质。网格性质包括:文本区域,图片区域,混合区域。混合区域是指网格内包含两种或两种以上不一致的内容,例如,包含有文本和图片,或者包含有背景不同的文本。如图3.2所示:

    3.2 用网格划分图像

    判断区域的规则如下:

    1)      文本区域的方差中等,但是跳变差比较大;

    2)      图片区域的跳变差比较小;

    3)      跳变差中等的区域,可以考虑为混合区域。一般地,混合区域方差比较大,而且其灰度均值与至少一个邻近区域差别比较大。

    4.          对于网格进行合并和拆分操作,各种性质非常接近的邻近网格,可以判断为内容一致的,应该进行合并。混合区域应进行拆分。

    如何进行拆分也是我们需要考虑的问题。实际应用中,合并操作可以省略,主要进行的其实是拆分操作。对于需要拆分的网格所对应的一小块图像,我们可以简单地利用边缘检测算子来寻找主要边缘。由于我们假定样张的不均匀背景是呈块状分布的,而且图片也是规则矩形,所以寻找到的边缘也应该具有直线或近似直线形状。这样就可以将原来的网格进一步切割成互不相交的、各部分基本一致的小矩形,如图3.3所示:


    3.3 对混合区域网格进行拆分

    3.3.2使用改进的条件熵方法进行阈值化

    对划分后的每个小区域使用条件熵方法进行阈值化,实际上是进行全局的阈值化。在这一步,我们要做的事情包括:

    1)    计算共存矩阵;

    2)    0,1,….L-1个灰度,利用共存矩阵分别计算以当前灰度值为阈值的总体条件熵值;

    3)    选择最大熵值所对应的灰度值为阈值。

    以下是计算共存矩阵、条件熵值的算法:

    算法一. 对给定的图像计算共存矩阵

    //算法名称:GetCooccurenceMatrix

    //输入:二维图像数据数组Image[n][n],灰度范围[0,L-1]

    //输出:共存矩阵二维数组Cmatrix[L][L]

    Begin:

        //将共存矩阵初始化,每个元素都置为0:

    Set each element in Cmatrix[L][L] as 0;

        //对图像中每个像素点,做下列操作:

        for(i = 0; i < n; i ++)

           for(j = 0; j < n; j ++)

           {

               //取当前像素点(i,j)的灰度值

               Set CurrentGrayLevel as Image[i][j];

               //遍历(i,j)的8个邻近像素点

               for(each Neighbour point in 8-neighbour of Image[i][j])

               {

                  //取其灰度值

                  Set NeighbourGrayLevel as Gray level of Neighbour point;

                  //将共存矩阵中对应位置的元素值加1:

    CMatrix[CurrentGrayLevel][ NeighbourGrayLevel] ++;

               }

           }

        return CMatrix[L][L];

    End;

     算法二. 利用共存矩阵求灰度值s对应的条件熵值:

    //算法名称:GetConditionalEntropy

    //输入:灰度值s,共存矩阵CMatrix[L][L]

    //输出:条件熵值H

    Begin:

        //初始化两个条件熵值H(O|B)和H(B|O)

    Set H(O|B) as 0;

    Set H(B|O) as 0;

        //计算 :

        for(i = 0; i < s; i ++)

           for(j = s+1; j <= L-1; j ++)

        {

        //首先计算:

           {

               temp = 0;

               for(l = 0; l < s; l ++)

           for(m = s+1; m <= L-1; m ++)

           {  

               temp += CMatrix[l][m];

           }

        po(i,j) = CMatrix[i][j]/temp;

           }

    //其次计算:

           {

               temp = 0;

               for(l = 0; l < s; l ++)

        {  

           temp += CMatrix[l][j];

        }

        po(i|j) = CMatrix[i][j]/temp;

           }

    //将计算得来的值加到H(O|B):

    H(O|B) += po(i,j) * exp(1 - po(i|j));

    }

    //计算 :

        for(i = s+1; i <= L-1; i ++)

           for(j = 0; j < s; j ++)

        {

        //首先计算:

           {

               temp = 0;

               for(l  = s+1; l<= L-1; l ++)

           for(m= 0; m< s; m ++)

           {  

               temp += CMatrix[l][m];

           }

        pb(i,j) = CMatrix[i][j]/temp;

           }

    //其次计算:

           {

               temp = 0;

               for(l = s+1; l <= L-1; l ++)

        {  

           temp += CMatrix[l][j];

        }

        pb(i|j) = CMatrix[i][j]/temp;

           }

    //将计算得来的值加到H(B|O):

    H(B|O) += pb(i,j) * exp(1 - pb(i|j));

    }  

    //总体条件熵值H为H(O|B) 与 H(B|O)之和:

    H = H(O|B) + H(B|O)

    return H;

    End;     

    算法三. 使用条件熵方法求图像分割阈值

    //算法名称:CoditionalEntropyMethod

    //输入:二维图像数据数组Image[n][n],灰度范围[0,L-1]

    //输出:最佳分割阈值T

    Begin:

        //调用算法GetCooccurenceMatrix,得到共存矩阵:

    CMatrix = GetCooccurenceMatrix (Image[n][n]);

        //对灰度值0,1,L-1,求最大条件熵值对应的灰度值:

        set MAX as 0;

        set T as 0;

        for(s = 0;s <= L-1; s ++)

        {

        //调用算法GetConditionalEntropy,求s对应的条件熵值H(s):

        H(s) = GetConditionalEntropy(s,CMatrix);

        if(H(s) > MAX)

        {

           MAX = H(s);

           T = s;

        }

    }

    return T;

    End;

    分析算法三,我们发现对它可以进行改进。本算法主要的时间空间消耗是在求最大熵值时对所有L个灰度值计算熵值并比较的循环中。而在实际图像中,阈值可能落在一个较小的灰度范围内。因此,我们可以在调用最大熵算法前,增加一个挑选可能的阈值集合的模块,然后在这个集合内应用最大熵选择操作,通过减少循环次数,可以大大减少计算量,而且不会对算法的效果造成影响。

    这个候选阈值集合应该如何挑选呢?我们注意到,如果区域主要由文本构成,那么其灰度分布直方图往往呈一大一小的双峰分布,其中较高的峰对应着背景,较低的峰对应着目标(文本)。记较高的峰对应的灰度值为,较低的峰对应的灰度值为,那么阈值显然应该落在之间,阈值寻找范围可以从[0,L-1]缩小到()。(不失一般性,可以假设<。)


    为了使直方图能够更好地体现灰度分布特性,我们在对其分析之前,可以先应用一个平滑操作,将原来比较粗糙的直方图变得平滑,消除毛刺、噪音、假“峰”等等的影响。

    3.4 a 原始灰度图像区域


    3.4 b 原始图像的灰度分布直方图

    3.4 c 经过平滑变换后的灰度分布直方图

    对平滑后的直方图寻找双峰,找出其对应的灰度值()。在此区间内使用最大条件熵算法寻找最佳阈值。

    3.4 d 在平滑后的直方图上寻找双峰,在其间计算阈值

    3.4 e 使用阈值分割后的二值图像

    改进后的阈值化算法流程如下:

    1)        获取各个小区域的灰度分布直方图;

    2)        对灰度直方图进行平滑,即用邻近点的频率均值代替当前频率值;

    3)        对于存在明显双峰的直方图对应的小区域,将阈值寻找范围缩小到双峰之间;

    4)        使用最大条件熵方法计算最佳阈值。

    3.4 算法效果比较

    实验结果表明,本文提出的算法对灰度变化范围比较大、背景呈块状变化的文档图像有比较好的二值化效果。图3.5显示了使用全局Otsu方法和用本算法进行阈值化的一些结果:     

    3.5 a原图像

    3.5 b 使用Otsu方法二值化的结果

    3.5 c 使用本文算法二值化的结果

    3.6 a 原图像

    3.6 b 使用Ostu方法二值化的结果

    3.6 c使用本文算法二值化的结果


    参考文献

    [1]     Pal,N.R.,Bhandari,D., “Image thresholding”. Signal Processing , 1993,33(2),pp.139-158

    [2]     C. K. Chow and T. Kaneko, "Automatic detection of the left ventricle from cineangiograms", Computers and Biomedical Research, , 1972,vol. 5, pp. 388—410

    [3]     J.Sauvola,M.Pietikäinen, “Adaptive document image binarization”, Pattern Recognition ,2000(33), pp. 225-236


    第四章 总结与展望

    在前面的讨论中,本文在对各种图像阈值化算法分析、研究的基础上提出了一种改进的灰值文档图像动态阈值化算法——直接局域二值化方法,从初步的实验结果来看,它对于背景与文字灰度不均匀的文档样张具有较好效果,比起原来的全局阈值方法有很大改进。

    但是,本算法还不是十分成熟,如果应用在某个实际系统中,在某些特殊情况下可能会出现一些问题。例如,在划分小区域时,如果某个小区域内只存在背景,且背景灰度有明显变化,可能在二值化以后就会出现明显的黑块,如第三章图3.6c中左方出现的黑块。但是一般来说,这样的错误对于文字区域并无影响,因此对以后的识别过程影响也不大,所以并不是致命的弱点。另外,当文档样张背景有纹理图案,特别是不规则纹理图案时,二值化后可能仍存在一些噪音点,这时可以考虑添加一个消除噪音模块。还有,本算法对于条件熵方法的优化处理,也还可以改进,使候选阈值集合进一步缩小。这里可以考虑借鉴第二章2.2.5.4中介绍的快速二维熵算法,对图像进行量化。这些问题都需要在以后的实际应用中,有针对性地对其加以改进。

    文档图像的阈值化是一个比较困难和复杂的问题,尤其是对于版面、颜色、纹理比较复杂的文档,怎样才能兼顾算法的效果与速度两方面,使之都在可忍受的范围内,是一个值得深入研究的领域。我们可以考虑引入MRF(马尔可夫随机场)模型来处理图像分割问题,以及使用模糊判别等等手段。这些都可以作为我们以后研究的方向

     



    [1]所谓Manhattan式版面是指文档图像中各对象可以用一个矩型区域来表示,这些区域互不重叠,而且之间具有明显的空白。

    转载于:https://www.cnblogs.com/kuaful/archive/2008/06/01/1211502.html

    展开全文
  • 图像拼接

    2015-11-24 22:07:17
    图像拼接的基本流程 (1) 图像预处理:对原始图像进行直方图匹配、平滑滤波、增强变换等数字图像 处理的基本操作,为图像拼接的下一步作好准备。 (2) 图像配准:图像配准是整个图像拼接流程的核心,...

    图像拼接的基本流程


    (1) 图像预处理:对原始图像进行直方图匹配、平滑滤波、增强变换等数字图像
    处理的基本操作,为图像拼接的下一步作好准备。

    (2) 图像配准:图像配准是整个图像拼接流程的核心,配准的精度决定了图像的拼接质量。其基本思想是:首先找到待配准图像与参考图像的模板或特征点的对应位置,然后根据对应关系建立参考图像与待配准图像之间的转换数学模型将待配准图像转换到参考图像的坐标系中,确定两图像之间的重叠区域。精确配准的关键是寻找一个能很好描述两幅图像转换关系的数据模型。

    (3) 图像合成:确定了两幅图像之间的转换关系模型,即重叠区域后,就需要根据重叠区域的信息将待拼接图像镶嵌成一个视觉可行的全景图。由于地形存在微小差别或拍摄条件不同等因素造成图像灰度(或亮度)差异,或者图像配准结果仍存在一定配准误差,为了尽可能地减少遗留变形或图像间的亮度(或灰度)差异对镶嵌结果的影响,就需要选择合适的图像合成策略。



    目前,国内外学者提出了很多图像配准的方法,但是各种方法都与一定范围的应用领域有关系,具有各自不同的特点。它们一般由四个要素组成[9]:特征空间、相似性度量、搜索空间和搜索策略。

    (1) 特征空间(feature space)
    特征空间是从图像中提取的用于配准的信息。特征可以是图像的灰度值,也可以是边界、轮廓等结构特征,或是角点、高曲率点等显著特征,亦或是统计特征、高层结构描述与句法描述等[9-11]。

    (2) 相似性度量(similarity metric)
    相似性度量是度量配准图像特征之间的相似性。典型的相似性度量有灰度相关、相关系数、互信息等。基于图像特征配准算法,常采用的相似性度量一般建立在各种距离函数上,如欧式距离、街区距离、Hausdorff 距离等。特征空间代表参与配准的数据,相似性度量决定配准度,二者的结合可以忽略许多与配准不相关的畸变,突出图像的本质结构与特征。

    (3) 搜索空间(search space)
    图像配准问题是一个参数的最优化估计问题,搜索空间就是指所有可能的变换组成的空间,即待估计参数组成的空间。搜索空间的组成和范围由图像变换的类型和强度决定。图像变换分为全局变换和局部变换。全局变换将整幅数字图像作为研究对象,用一个参数矩阵来描述整个图像的变换参数,常见的全局几何变换有仿射变换、投影变换和非线性变换等;局部变换允许变换参数有位置依赖性,即各个单元的变换参数随其所处位置不同而不同。配准算法就是要在搜索空间找出一个使图像之间的相似性度量最佳的位置。

    (4) 搜索策略(search strategy)
    搜索策略是指在搜索空间中采用合适的方法找出平移、旋转等变换参数的最优估计。这对于减少配准特征和相似性度量的计算量具有重要意义。图像间的变换越复杂,搜索空间就越复杂,对搜索策略的要求就越高,因此选择合适的探索策略至关重要。常见的搜索策略有穷尽搜索、启发式搜索、广义Hough 变换、多尺度搜索、树与图搜索、序贯判决、线性规划、神经网络、遗传算法、模拟退火算法等等。每种方法都有它的优缺点,在很大程度上,搜索策略的选择取决于搜索空间的特性。

    设计配准算法时首先确定图像的类型和成像畸变范围,然后根据图像配准的性能指标确定特征空间和搜索空间,最后选择合适的搜索策略找到能使待配准图像之间的相似性度量最大的最佳匹配关系模型。根据上述四个要素的不同选择,产生了各种具体的图像配准技术的不同分类方法。目前提出的用于图像配准的算法可以分为三大类型:基于区域灰度相关的匹配算法、基于特征相关的匹配算法和基于解释相似的算法[12]。图像配准技术经过多年的研究,每类方法中都包含了不同的具体实现途径以适应具体问题。其中,基于解释的图像匹配需要建立在图片自动判读的专家系统之上进行,至今尚未取得突破性的进展。

    —————————————————————————————————————

    作者:刘伟民

    毕业于:中国科学院计算技术研究所   

     

     

       最近我在看图像拼接的东西,在这里普及一下知识。

        图像拼接就是将两幅图像(有公共部分)拼成一个大的图像。这个用在视频里面就是视频拼接。视频拼接中有着特殊的需求:

    1. 视频中的图像一般重叠面积比较大,因此容易找到对应点。
    2. 视频中的图像一般很难用到和人交互的东西,所以一般识别对应点应该是自动的过程。
    3. 视频拼接一般要求速度比较快,最好达到实时拼接的效果。

    图像拼接的要求:

    1. 很多图像要自动识别出图像之间的序列关系,就是那幅图像应该和那副图像挨着。
    2. 图像不要求速度,但是要求精确,就是准确度比较高。
    3. 图像拼接可以加入人工成分,就是可以让人自己选定对应的点。

        在图像选完了对应的特征点之后就该进行图像对应特征点的计算,求出两幅图像的对应的透视变换关系。

        变换

        接下来进行拼接,将两幅图像合成一幅图像。

     

     

    下面我们用SIFT特征进行一个简单的图像拼接过程:

    步骤:

     提取图像特征

     对特征进行匹配,找到正确的匹配点

     根据这些特征对,进行透视变换矩阵计算。

     根据得到的透视变换矩阵,对两幅图像进行透视映射。

     将变换后的两幅图像加到一个图像中。

    拼接实验:

    素材图片:

     

    基于sift特征的简易图像拼接程序

    基于sift特征的简易图像拼接程序

    基于sift特征的简易图像拼接程序

    基于sift特征的简易图像拼接程序

     

    下面是我拼接的结果:



    基于sift特征的简易图像拼接程序

     

    下面是*****的拼接结果:



    基于sift特征的简易图像拼接程序

    结果分析:

    我的实验中的图像有的地方很亮,这个对图像进行加权融合就可以解决。

    其效果和*****的差不多。

    *****的拼接程序经常出现【拼接失败】的情况,我的程序更具有鲁棒性。

    缺点:拼接速度没有*****的快。

     

    基于sift特征的简易图像拼接程序

    基于sift特征的简易图像拼接程序

    基于sift特征的简易图像拼接程序

    基于sift特征的简易图像拼接程序

    基于sift特征的简易图像拼接程序 

    基于sift特征的简易图像拼接程序

    ——————————————————————————————————————
    算法描述 
    
    procedure ImageMatching 
    
    { 
    
        输入FirstImage; 
    
        输入SecondImage; 
    
        //获得两幅图象的大小 
    
        Height1=GetImageHeight(FirstImage); 
    
        Height2=GetImageHeight(SecondImage); 
    
        Width1=GetImageWidth(FirstImage); 
    
        Width2=GetImageWidth(SecondImage); 
    
    // 从第二幅图象取网格匹配模板 
    
        SecondImageGrid = GetSecondImageGrid(SecondImage); 
    
    // 粗略匹配,网格在第一幅图象中先从左向右移动,再从下到上移动,每次移动一个网格间距,Step_Width 或Step_Height,当网格移出重叠区域后结束 
    
    y=Heitht1-GridHeight; 
    
    MinValue = MaxInteger; 
    
    While ( y<Height1-OverlapNumber)//当网格移出重叠部分后结束 
    
    { 
    
           x=Grid_Width/2; //当网格位于第一幅图象的最左边时,A点的横坐标。 
    
           While ( x<(Width1-Grid_Width/2) ) 
    
    { 
    
               FirstImageGrid=GetImgaeGrid(FirstImgaeGrid, x, y); 
    
            differ=CaculateDiff(FirstImgaeGrid, SecondImageGrid);//计算象素值差的平 
    
                                   //方和 
    
               if (differ<MinValue) 
    
               { 
    
                  BestMatch_x=x; 
    
                  BestMatch_y=y; 
    
                  MinValue = differ; 
    
               } 
    
           x= x+Step_width; 
    
           } 
    
           y=y-Step_Height; 
    
    }   
    
    //精确匹配 
    
    Step_Width= Step_Width/2; 
    
    Step_Height= Step_Height/2; 
    
    While ( Step_Height>0 & Step_Width>0)//当水平步长和垂直步长均减为零时结束 
    
    { 
    
        if(Step_Height==0)//当仅有垂直步长减为零时,将其置为1 
    
           Step_Height=1; 
    
        If(Step_Width==0) //当仅有水平步长减为零时,将其置为1 
    
           Step_Width=1; 
    
    temp_x = BestMatch_x; 
    
    temp_y = BestMatch_y; 
    
        for ( i= -1; i<1; i++) 
    
           for( j= -1; j<1; j++) 
    
           { 
    
               if ((i=0&j!=0)|(i!=0&j=0)) 
    
               { 
    
                  FirstImageGrid=GetImgaeGrid(FirstImgaeGrid, 
    
    temp_x+i*Step_Width, temp_y +j*Step_Height); 
    
               differ=CaculateDiff(FirstImgaeGrid, SecondImageGrid); 
    
           if (differ<MinValue) 
    
                  { 
    
                      BestMatch_x=x; 
    
                      BestMatch_y=y; 
    
                      MinValue = differ; 
    
               } 
    
               } 
    
              } 
    
           Step_Height = Step_Height /2; 
    
           Step_Width = Step_Width/2; 
    
        } 
    
    }  

    ——————————————————————————————————————

    主要分为以下几个步骤:

    (1) 读入两张图片并分别提取SIFT特征

    (2) 利用k-d tree和BBF算法进行特征匹配查找

    (3) 利用RANSAC算法筛选匹配点并计算变换矩阵

    (3) 图像融合


    SIFT算法以及RANSAC算法都是利用的RobHess的SIFT源码,前三个步骤RobHess的源码中都有自带的示例。


    (1) SIFT特征提取

    直接调用RobHess源码(RobHess的SIFT源码分析:综述)中的sift_features()函数进行默认参数的SIFT特征提取,主要代码如下:

    1. img1_Feat = cvCloneImage(img1);//复制图1,深拷贝,用来画特征点  
    2. img2_Feat = cvCloneImage(img2);//复制图2,深拷贝,用来画特征点  
    3.   
    4. //默认提取的是LOWE格式的SIFT特征点  
    5. //提取并显示第1幅图片上的特征点  
    6. n1 = sift_features( img1, &feat1 );//检测图1中的SIFT特征点,n1是图1的特征点个数  
    7. export_features("feature1.txt",feat1,n1);//将特征向量数据写入到文件  
    8. draw_features( img1_Feat, feat1, n1 );//画出特征点  
    9. cvNamedWindow(IMG1_FEAT);//创建窗口  
    10. cvShowImage(IMG1_FEAT,img1_Feat);//显示  
    11.   
    12. //提取并显示第2幅图片上的特征点  
    13. n2 = sift_features( img2, &feat2 );//检测图2中的SIFT特征点,n2是图2的特征点个数  
    14. export_features("feature2.txt",feat2,n2);//将特征向量数据写入到文件  
    15. draw_features( img2_Feat, feat2, n2 );//画出特征点  
    16. cvNamedWindow(IMG2_FEAT);//创建窗口  
    17. cvShowImage(IMG2_FEAT,img2_Feat);//显示  
    检测出的SIFT特征点如下:

                     


    (2) 利用k-d tree和BBF算法进行特征匹配查找,并根据最近邻和次近邻距离比值进行初步筛选

    也是调用RobHess源码中的函数,加上之后的一些筛选处理,主要代码如下:

    1. //根据图1的特征点集feat1建立k-d树,返回k-d树根给kd_root  
    2. kd_root = kdtree_build( feat1, n1 );  
    3.   
    4. Point pt1,pt2;//连线的两个端点  
    5. double d0,d1;//feat2中每个特征点到最近邻和次近邻的距离  
    6. int matchNum = 0;//经距离比值法筛选后的匹配点对的个数  
    7.   
    8. //遍历特征点集feat2,针对feat2中每个特征点feat,选取符合距离比值条件的匹配点,放到feat的fwd_match域中  
    9. for(int i = 0; i < n2; i++ )  
    10. {  
    11.     feat = feat2+i;//第i个特征点的指针  
    12.     //在kd_root中搜索目标点feat的2个最近邻点,存放在nbrs中,返回实际找到的近邻点个数  
    13.     int k = kdtree_bbf_knn( kd_root, feat, 2, &nbrs, KDTREE_BBF_MAX_NN_CHKS );  
    14.     if( k == 2 )  
    15.     {  
    16.         d0 = descr_dist_sq( feat, nbrs[0] );//feat与最近邻点的距离的平方  
    17.         d1 = descr_dist_sq( feat, nbrs[1] );//feat与次近邻点的距离的平方  
    18.         //若d0和d1的比值小于阈值NN_SQ_DIST_RATIO_THR,则接受此匹配,否则剔除  
    19.         if( d0 < d1 * NN_SQ_DIST_RATIO_THR )  
    20.         {   //将目标点feat和最近邻点作为匹配点对  
    21.             pt2 = Point( cvRound( feat->x ), cvRound( feat->y ) );//图2中点的坐标  
    22.             pt1 = Point( cvRound( nbrs[0]->x ), cvRound( nbrs[0]->y ) );//图1中点的坐标(feat的最近邻点)  
    23.             pt2.x += img1->width;//由于两幅图是左右排列的,pt2的横坐标加上图1的宽度,作为连线的终点  
    24.             cvLine( stacked, pt1, pt2, CV_RGB(255,0,255), 1, 8, 0 );//画出连线  
    25.             matchNum++;//统计匹配点对的个数  
    26.             feat2[i].fwd_match = nbrs[0];//使点feat的fwd_match域指向其对应的匹配点  
    27.         }  
    28.     }  
    29.     free( nbrs );//释放近邻数组  
    30. }  
    31. //显示并保存经距离比值法筛选后的匹配图  
    32. cvNamedWindow(IMG_MATCH1);//创建窗口  
    33. cvShowImage(IMG_MATCH1,stacked);//显示  
    匹配结果如下:



    (3) 利用RANSAC算法筛选匹配点并计算变换矩阵

    此部分最主要的是RobHess源码中的ransac_xform()函数,此函数实现了用RANSAC算法筛选匹配点,返回结果是计算好的变换矩阵。

    此部分中,我利用匹配点的坐标关系,对输入的两幅图像的左右关系进行了判断,并根据结果选择使用矩阵H或H的逆阵进行变换。

    所以读入的两幅要拼接的图像的左右位置关系可以随意,程序中可自动调整。

    主要代码如下:


    1. //利用RANSAC算法筛选匹配点,计算变换矩阵H,  
    2. //无论img1和img2的左右顺序,计算出的H永远是将feat2中的特征点变换为其匹配点,即将img2中的点变换为img1中的对应点  
    3. H = ransac_xform(feat2,n2,FEATURE_FWD_MATCH,lsq_homog,4,0.01,homog_xfer_err,3.0,&inliers,&n_inliers);  
    4.   
    5. //若能成功计算出变换矩阵,即两幅图中有共同区域  
    6. if( H )  
    7. {  
    8.     qDebug()<<tr("经RANSAC算法筛选后的匹配点对个数:")<<n_inliers<<endl; //输出筛选后的匹配点对个数  
    9.   
    10.     int invertNum = 0;//统计pt2.x > pt1.x的匹配点对的个数,来判断img1中是否右图  
    11.   
    12.     //遍历经RANSAC算法筛选后的特征点集合inliers,找到每个特征点的匹配点,画出连线  
    13.     for(int i=0; i<n_inliers; i++)  
    14.     {  
    15.         feat = inliers[i];//第i个特征点  
    16.         pt2 = Point(cvRound(feat->x), cvRound(feat->y));//图2中点的坐标  
    17.         pt1 = Point(cvRound(feat->fwd_match->x), cvRound(feat->fwd_match->y));//图1中点的坐标(feat的匹配点)  
    18.         //qDebug()<<"pt2:("<<pt2.x<<","<<pt2.y<<")--->pt1:("<<pt1.x<<","<<pt1.y<<")";//输出对应点对  
    19.   
    20.         //统计匹配点的左右位置关系,来判断图1和图2的左右位置关系  
    21.         if(pt2.x > pt1.x)  
    22.             invertNum++;  
    23.   
    24.         pt2.x += img1->width;//由于两幅图是左右排列的,pt2的横坐标加上图1的宽度,作为连线的终点  
    25.         cvLine(stacked_ransac,pt1,pt2,CV_RGB(255,0,255),1,8,0);//在匹配图上画出连线  
    26.     }  
    27.   
    28.     cvNamedWindow(IMG_MATCH2);//创建窗口  
    29.     cvShowImage(IMG_MATCH2,stacked_ransac);//显示经RANSAC算法筛选后的匹配图  
    30.   
    31.     /*程序中计算出的变换矩阵H用来将img2中的点变换为img1中的点,正常情况下img1应该是左图,img2应该是右图。 
    32.       此时img2中的点pt2和img1中的对应点pt1的x坐标的关系基本都是:pt2.x < pt1.x 
    33.       若用户打开的img1是右图,img2是左图,则img2中的点pt2和img1中的对应点pt1的x坐标的关系基本都是:pt2.x > pt1.x 
    34.       所以通过统计对应点变换前后x坐标大小关系,可以知道img1是不是右图。 
    35.       如果img1是右图,将img1中的匹配点经H的逆阵H_IVT变换后可得到img2中的匹配点*/  
    36.   
    37.     //若pt2.x > pt1.x的点的个数大于内点个数的80%,则认定img1中是右图  
    38.     if(invertNum > n_inliers * 0.8)  
    39.     {  
    40.         CvMat * H_IVT = cvCreateMat(3, 3, CV_64FC1);//变换矩阵的逆矩阵  
    41.         //求H的逆阵H_IVT时,若成功求出,返回非零值  
    42.         if( cvInvert(H,H_IVT) )  
    43.         {  
    44.             cvReleaseMat(&H);//释放变换矩阵H,因为用不到了  
    45.             H = cvCloneMat(H_IVT);//将H的逆阵H_IVT中的数据拷贝到H中  
    46.             cvReleaseMat(&H_IVT);//释放逆阵H_IVT  
    47.             //将img1和img2对调  
    48.             IplImage * temp = img2;  
    49.             img2 = img1;  
    50.             img1 = temp;  
    51.             ui->mosaicButton->setEnabled(true);//激活全景拼接按钮  
    52.         }  
    53.         else//H不可逆时,返回0  
    54.         {  
    55.             cvReleaseMat(&H_IVT);//释放逆阵H_IVT  
    56.             QMessageBox::warning(this,tr("警告"),tr("变换矩阵H不可逆"));  
    57.         }  
    58.     }  
    59.     else  
    60.         ui->mosaicButton->setEnabled(true);//激活全景拼接按钮  
    61. }  
    62. else //无法计算出变换矩阵,即两幅图中没有重合区域  
    63. {  
    64.     QMessageBox::warning(this,tr("警告"),tr("两图中无公共区域"));  
    65. }  
    经RANSAC筛选后的匹配结果如下图:



    (3) 图像融合

    这里有两种拼接方法:

    ① 简易拼接方法的过程是:首先将右图img2经变换矩阵H变换到一个新图像中,然后直接将左图img1加到新图像中,这样拼接出来会有明显的拼接缝,但也是一个初步的成品了。

    ② 另一种方法首先也是将右图img2经变换矩阵H变换到一个新图像中,然后图像的融合过程将目标图像分为三部分,最左边完全取自img1中的数据,中间的重合部分是两幅图像的加权平均,重合区域右边的部分完全取自img2经变换后的图像。加权平均的权重选择也有好多方法,比如可以使用最基本的取两张图像的平均值,但这样会有明显的拼接缝。这里首先计算出拼接区域的宽度,设d1,d2分别是重叠区域中的点到重叠区域左边界和右边界的距离,则使用如下公式计算重叠区域的像素值:

    ,这样就可以实现平滑过渡。

    主要代码如下:

    1. //若能成功计算出变换矩阵,即两幅图中有共同区域,才可以进行全景拼接  
    2. if(H)  
    3. {  
    4.     //拼接图像,img1是左图,img2是右图  
    5.     CalcFourCorner();//计算图2的四个角经变换后的坐标  
    6.     //为拼接结果图xformed分配空间,高度为图1图2高度的较小者,根据图2右上角和右下角变换后的点的位置决定拼接图的宽度  
    7.     xformed = cvCreateImage(cvSize(MIN(rightTop.x,rightBottom.x),MIN(img1->height,img2->height)),IPL_DEPTH_8U,3);  
    8.     //用变换矩阵H对右图img2做投影变换(变换后会有坐标右移),结果放到xformed中  
    9.     cvWarpPerspective(img2,xformed,H,CV_INTER_LINEAR + CV_WARP_FILL_OUTLIERS,cvScalarAll(0));  
    10.     cvNamedWindow(IMG_MOSAIC_TEMP); //显示临时图,即只将图2变换后的图  
    11.     cvShowImage(IMG_MOSAIC_TEMP,xformed);  
    12.   
    13.     //简易拼接法:直接将将左图img1叠加到xformed的左边  
    14.     xformed_simple = cvCloneImage(xformed);//简易拼接图,克隆自xformed  
    15.     cvSetImageROI(xformed_simple,cvRect(0,0,img1->width,img1->height));  
    16.     cvAddWeighted(img1,1,xformed_simple,0,0,xformed_simple);  
    17.     cvResetImageROI(xformed_simple);  
    18.     cvNamedWindow(IMG_MOSAIC_SIMPLE);//创建窗口  
    19.     cvShowImage(IMG_MOSAIC_SIMPLE,xformed_simple);//显示简易拼接图  
    20.   
    21.     //处理后的拼接图,克隆自xformed  
    22.     xformed_proc = cvCloneImage(xformed);  
    23.   
    24.     //重叠区域左边的部分完全取自图1  
    25.     cvSetImageROI(img1,cvRect(0,0,MIN(leftTop.x,leftBottom.x),xformed_proc->height));  
    26.     cvSetImageROI(xformed,cvRect(0,0,MIN(leftTop.x,leftBottom.x),xformed_proc->height));  
    27.     cvSetImageROI(xformed_proc,cvRect(0,0,MIN(leftTop.x,leftBottom.x),xformed_proc->height));  
    28.     cvAddWeighted(img1,1,xformed,0,0,xformed_proc);  
    29.     cvResetImageROI(img1);  
    30.     cvResetImageROI(xformed);  
    31.     cvResetImageROI(xformed_proc);  
    32.     cvNamedWindow(IMG_MOSAIC_BEFORE_FUSION);  
    33.     cvShowImage(IMG_MOSAIC_BEFORE_FUSION,xformed_proc);//显示融合之前的拼接图  
    34.   
    35.     //采用加权平均的方法融合重叠区域  
    36.     int start = MIN(leftTop.x,leftBottom.x) ;//开始位置,即重叠区域的左边界  
    37.     double processWidth = img1->width - start;//重叠区域的宽度  
    38.     double alpha = 1;//img1中像素的权重  
    39.     for(int i=0; i<xformed_proc->height; i++)//遍历行  
    40.     {  
    41.         const uchar * pixel_img1 = ((uchar *)(img1->imageData + img1->widthStep * i));//img1中第i行数据的指针  
    42.         const uchar * pixel_xformed = ((uchar *)(xformed->imageData + xformed->widthStep * i));//xformed中第i行数据的指针  
    43.         uchar * pixel_xformed_proc = ((uchar *)(xformed_proc->imageData + xformed_proc->widthStep * i));//xformed_proc中第i行数据的指针  
    44.         for(int j=start; j<img1->width; j++)//遍历重叠区域的列  
    45.         {  
    46.             //如果遇到图像xformed中无像素的黑点,则完全拷贝图1中的数据  
    47.             if(pixel_xformed[j*3] < 50 && pixel_xformed[j*3+1] < 50 && pixel_xformed[j*3+2] < 50 )  
    48.             {  
    49.                 alpha = 1;  
    50.             }  
    51.             else  
    52.             {   //img1中像素的权重,与当前处理点距重叠区域左边界的距离成正比  
    53.                 alpha = (processWidth-(j-start)) / processWidth ;  
    54.             }  
    55.             pixel_xformed_proc[j*3] = pixel_img1[j*3] * alpha + pixel_xformed[j*3] * (1-alpha);//B通道  
    56.             pixel_xformed_proc[j*3+1] = pixel_img1[j*3+1] * alpha + pixel_xformed[j*3+1] * (1-alpha);//G通道  
    57.             pixel_xformed_proc[j*3+2] = pixel_img1[j*3+2] * alpha + pixel_xformed[j*3+2] * (1-alpha);//R通道  
    58.         }  
    59.     }  
    60.     cvNamedWindow(IMG_MOSAIC_PROC);//创建窗口  
    61.     cvShowImage(IMG_MOSAIC_PROC,xformed_proc);//显示处理后的拼接图  
    62.   
    63.     /*重叠区域取两幅图像的平均值,效果不好 
    64.         //设置ROI,是包含重叠区域的矩形 
    65.         cvSetImageROI(xformed_proc,cvRect(MIN(leftTop.x,leftBottom.x),0,img1->width-MIN(leftTop.x,leftBottom.x),xformed_proc->height)); 
    66.         cvSetImageROI(img1,cvRect(MIN(leftTop.x,leftBottom.x),0,img1->width-MIN(leftTop.x,leftBottom.x),xformed_proc->height)); 
    67.         cvSetImageROI(xformed,cvRect(MIN(leftTop.x,leftBottom.x),0,img1->width-MIN(leftTop.x,leftBottom.x),xformed_proc->height)); 
    68.         cvAddWeighted(img1,0.5,xformed,0.5,0,xformed_proc); 
    69.         cvResetImageROI(xformed_proc); 
    70.         cvResetImageROI(img1); 
    71.         cvResetImageROI(xformed); //*/  
    72.   
    73.     /*对拼接缝周围区域进行滤波来消除拼接缝,效果不好 
    74.         //在处理前后的图上分别设置横跨拼接缝的矩形ROI 
    75.         cvSetImageROI(xformed_proc,cvRect(img1->width-10,0,img1->width+10,xformed->height)); 
    76.         cvSetImageROI(xformed,cvRect(img1->width-10,0,img1->width+10,xformed->height)); 
    77.         cvSmooth(xformed,xformed_proc,CV_MEDIAN,5);//对拼接缝周围区域进行中值滤波 
    78.         cvResetImageROI(xformed); 
    79.         cvResetImageROI(xformed_proc); 
    80.         cvShowImage(IMG_MOSAIC_PROC,xformed_proc);//显示处理后的拼接图 */  
    81.   
    82.     /*想通过锐化解决变换后的图像失真的问题,对于扭曲过大的图像,效果不好 
    83.         double a[]={  0, -1,  0, -1,  5, -1, 0, -1,  0  };//拉普拉斯滤波核的数据 
    84.         CvMat kernel = cvMat(3,3,CV_64FC1,a);//拉普拉斯滤波核 
    85.         cvFilter2D(xformed_proc,xformed_proc,&kernel);//滤波 
    86.         cvShowImage(IMG_MOSAIC_PROC,xformed_proc);//显示处理后的拼接图*/  
    87.   
    88. }  

    右图经变换后的结果如下图:


    简易拼接结果如下图:


    使用第二种方法时,重合区域融合之前如下图:


    加权平均融合之后如下图:



    用Qt做了个简单的界面,如下:



    还有很多不足,经常有黑边无法去除,望大家多多指正。

    实验环境:

    IDE是 Qt Creator 2.4.1Qt版本是4.8.4OpenCV版本是2.4.4,SIFT源码来自RobHess给出的源码

    源码下载:基于SIFT特征的全景图像拼接,Qt工程:http://download.csdn.net/detail/masikkk/5702681

    ——————————————————————————————————————

    遥感图像拼接与裁剪

    一、试验目的:通过ENVI和ARCGIS软件对下载的大足地区的TM图像进行拼接和裁剪,将自己的家乡所在的区域裁剪出来。通过本次试验,初步熟悉ENVI和ARCGIS软件,为今后环境遥感学习奠定基础。

    二、试验原理:  1、图像拼接原理:ENVI的图像拼接功能提供交互式的方式将没有地理坐标或者有地理坐标的多幅图像合并,生成一幅单一的合成图像。2、图像裁剪原理:通常按照行政区划边界或自然区划边界进行图像剪裁,在基础数据生产中,还经常要进行标准分幅裁剪。

    三、数据来源:已经数字化好的1:400万部全国各县界的shp格式的矢量数据。从遥感数据共享中心(http://ids.ceode.ac.cn/)下载下来的LANDSAT5 TM数据。波段数为7,各波段信息如表1所示:

    波段

    波长范围(μm)

    分辨率

    1

    0.45~0.53

    30米

    2

    0.52~0.60

    30米

    3

    0.63~0.69

    30米

    4

    0.76~0.90

    30米

    5

    1.55~1.75

    30米

    6

    10.40~12.50

    120>米

    7

    2.08~2.35

    30米

    表1

    四、实验步骤:

    1、数据的下载:

        打开遥感数据共享中心(http://ids.ceode.ac.cn/)并登陆用户。

    [转载]遥感图像预处理之图像拼接与裁剪
    图1

        按空间条件进行查询,输入中国、重庆市、重庆市郊、大足县等信息,点击显示查询所需下载的区域。并且选择LANDSAT5卫星进行查询。

    [转载]遥感图像预处理之图像拼接与裁剪

     

    图2

        通过查看JP2K快视图可以看出,并没有一副遥感影像能够完全覆盖所选的大足县,因此需要下载两幅图像进行图像的拼接。

    [转载]遥感图像预处理之图像拼接与裁剪
    图3

    2、图像的拼接

        打开ENVI软件,并通过3-2-1真彩色通道打开下载好的两幅遥感影像。

    [转载]遥感图像预处理之图像拼接与裁剪

    图4

        在ENVI主菜单中选择Map—Mosaicing—Georeferenced,在Mosaic对话框中点击Import—Import File,选择需要拼接的两幅图像。

    [转载]遥感图像预处理之图像拼接与裁剪
    图5

    [转载]遥感图像预处理之图像拼接与裁剪

     

    图6

        对图片点击右键,选中Edit Entry,在Edit Entry对话框中,设置Data Value to Ignor:0,忽略0值,设置Feathering Distance为10,羽化半径为10个像素,点击OK确定。

    [转载]遥感图像预处理之图像拼接与裁剪

     

    图8

        如图8所示,图像已经无接缝地拼接成功,点击File—Apply输出图片。

        将保存好的图片转存成TIFF文件,以便可以在ARCMAP中存在。在ENVI主菜单点击FILE—SAVE FILE AS—TIFF/GeoTIFF。

    [转载]遥感图像预处理之图像拼接与裁剪

     

    图9

    3、投影转换

        打开ARCMAP,导入中国县界数据和已经拼接好的遥感数据,通过字段查找(Find功能)查询属性字段为大足县的面数据。右键—Select选中大足县。

    [转载]遥感图像预处理之图像拼接与裁剪

    图10

        在ArcToolbox中,选择分析工具中的提取工具中的Clip裁剪,将所需的大足县面数据裁剪出来。

    [转载]遥感图像预处理之图像拼接与裁剪

    图11

    [转载]遥感图像预处理之图像拼接与裁剪
    图12

        同时打开拼接好的遥感影像和大足矢量文件可以发现,根据地理坐标,两幅图叠加在一起,矢量文件的位置就是我们需要裁剪的位置。

        在TOC工具栏里,对大足矢量文件点击右键—properties,并点击source选项卡查看其属性,发现其为兰伯特投影。如图13所示。

    [转载]遥感图像预处理之图像拼接与裁剪
    图13

        而打开拼接好的遥感影像却发现,其为横轴麦卡托投影。如图14所示:

    [转载]遥感图像预处理之图像拼接与裁剪

    图14

        由于ENVI软件和ARCGIS软件不同,ENVI软件中没有自动投影转换功能,因此,我们需要对这两幅图像进行投影转换,否则这两幅图像在ENVI软件中无法进行叠加。就无法进一步对图像进行裁剪。

        麦卡托投影是低纬度地区常用的一种投影方式。它将地球表面上各点投影到与地球相切的圆柱面上,然后依某条母线剪开而成。它是一种保角圆柱面投影。其计算公式如下:

    [转载]遥感图像预处理之图像拼接与裁剪

     

        兰勃特投影以圆锥为投影面构成保角投影,此圆锥与地球上某一纬度线相切或与两条纬度线相割,投影后沿某一母线剖开展平即成为兰勃特投影。

     

    [转载]遥感图像预处理之图像拼接与裁剪

     

        将兰伯特投影转换成麦卡托投影就是要将兰伯特投影公式中的X和Y倒推成地理坐标,然后再代入麦卡托投影公式求出麦卡托投影中新的X和Y。

        这些计算比较复杂,在ArcMap中,用封装好的工具很容易将其实现。在ArctoolBOX中,点击Data Management Tools—Projections and Transformations—Feature—Project。

    [转载]遥感图像预处理之图像拼接与裁剪

     

    图15

        我们需要将大足县矢量数据的兰伯特投影转换为横轴麦卡托投影(与大足裁剪一致),因此在Project对话框中,输入大足县矢量数据,在output dataset or feature class,点击右边的输入按钮,在Spatial Reference Properties点击Import选择拼接好的大足遥感影像数据。这部操作的含义是,将大足矢量数据的投影转换成为大足遥感影像数据的投影,并不需要我们知道转换过程,点击确定。

    [转载]遥感图像预处理之图像拼接与裁剪

    图16

        如图所示,转换后的大足县shp文件已经变成了横轴麦卡托投影,与遥感图像的投影一致了。

    4、图像的裁剪

        ENVI软件中图像的裁剪有三种方法:第一种为规则裁剪,第二种为选区ROI感兴趣区域进行掩膜裁剪,第三种就是我们首先要介绍的矢量文件掩膜裁剪。

    (1)矢量文件掩膜裁剪

        打开ENVI软件,输入拼接后的遥感图像。在ENVI主菜单中点击FILE—Open Vector File输入矢量文件(投影转换后的大足shp文件)。

    [转载]遥感图像预处理之图像拼接与裁剪

    图17

        在Available Vector List对话框中,选择FILE—Export layers to ROIs,点击大足接好图片,选择Convert all records of an EVF layer to one ROI,回到Available Vector List对话框,点击Load Selected可以查看ROI区域是否在图像指定位置。点击Load Selected之后,选择显示的窗口为遥感图像打开的窗口,将Current Layer颜色改成黑色。点击Apply查看,效果如图18所示:

    [转载]遥感图像预处理之图像拼接与裁剪

     

    图18

        当确定所需裁剪的区域是黑色框体选中的位置之后,在ENVI主菜单选择Basic Tools—Subset data via ROIs,选择大足拼接图像,在Spatial subset via ROI Parameters对话框中,选择刚才导出的ROI,将掩膜项设为YES,掩膜背景的值为0,点击确定输出文件。

    裁剪好的大足图像如图19所示:

    [转载]遥感图像预处理之图像拼接与裁剪

    图19

    (2)规则裁剪

        规则裁剪裁剪出来的区域一般为矩形,规则裁剪一般用于已经明确知道需要裁剪的区域范围所在,或者只需要大概裁剪出部分矩形区域,以缩小处理的面积。具体的ENVI操作步骤如下。

        在ENVI主菜单中点击FILE—SAVE FILE AS—ENVI STANDARD,如图20所示。在NEW FILE BUILDER对话框中,点击IMPORT FILE按钮,在CREATE NEW FILE INPUT FILE对话框中,选择大足接好,点击下方的spatial subset按钮。在Select Spatial Subset对话框中,可以通过选择裁剪的行和列来裁剪图像,如图21所示。也可以点击Image,通过点击红色方框的四个角来扩大缩小裁剪区域,通过点击方框内部实现裁剪区域的移动,也可以实现裁剪区域的选取,如图22所示。

    [转载]遥感图像预处理之图像拼接与裁剪
    图20

    [转载]遥感图像预处理之图像拼接与裁剪

     

    图21

    [转载]遥感图像预处理之图像拼接与裁剪

     

    图22

        点击OK输出图像即为所裁剪的矩形规则区域。如如图23所示:

    [转载]遥感图像预处理之图像拼接与裁剪

     

    图23

     

    (3)感兴趣区掩膜裁剪

        这种裁剪方法常常用于各类地物的提取裁剪,如水体等。

        选取感兴趣区的方法多种多样,如监督分类、非监督分类、NDVI和自定义波段运算、手工勾勒等等。

        下面就以大足拼接遥感图像中某区域的森林提取为例子,来说明这种感兴趣区掩膜裁剪的方法。

        首先对图像进行NDVI计算。严格意义上定量计算NDVI植被指数的值,需要将原始图像(DN值)经过辐射定标、大气矫正后得到地表反射率才能进行NDVI计算。但是在此为了简单起见,在不需要定量计算植被指数的时候,直接用DN值进行NDVI计算。

        如图24所示,在ENVI主菜单中,选择Transform—NDVI,我们用的是TM数据,对应的近红外波段为4波段,红光波段为3波段,ENVI软件中已经集成好,不需要我们进行设置。

    [转载]遥感图像预处理之图像拼接与裁剪

    图24

        NDVI计算结果如图25所示,在图中可以发现,植被区域,NDVI的值比较大,约为0.65—0.8,我们可以对此设置波段阈值。在ENVI主菜单中,选择Basic Tools—Region of Interest—Band threshold to ROI,选择刚才计算的NDVI,设置阈值为0.65—0.8如图27所示。

    [转载]遥感图像预处理之图像拼接与裁剪

    图25

    [转载]遥感图像预处理之图像拼接与裁剪

    图26

    [转载]遥感图像预处理之图像拼接与裁剪

     

    图27

        在ENVI主菜单选择Basic Tools—Subset data via ROIs,选择大足拼接图像,在Spatial subset via ROI Parameters对话框中,选择刚才NDVI提取的植被ROI,然后进行掩膜裁剪,其操作方法和矢量文件掩膜裁剪一致,在此不再赘述。

    五、实验结果与分析:

        经上诉步骤,通过不同的方法可以对图像进行不同的裁剪,最终获得了较好的效果。

    六、实验心得:

        实验操作需要对软件通过不断的练习与摸索才能比较熟练,清楚地了解到软件操作所带来的对自己有价值的东西。

    ——————————————————————————————————————
     经过数字化扫描及几何校正后的数字化遥感影像,均为一幅幅具有相同比例尺的影像图。这些影像图互相之间都存在着部分的重叠。所谓图像拼接就是通过对相邻影像图的无缝拼接处理,把这些影像图相互间的重叠部分去掉,从而为在逻辑上将这些影像图整合成覆盖区域的一幅影像图创造条件。

        图像拼接的具体工作步骤为: 

        首先是进行色差处理,借助photoshop软件中的色彩调整功能,将需要拼接的两幅相邻影像图的色彩调整到尽可能和谐。其次是选择拼接线,在两幅相邻的影像图上,用彩色线把需要进行拼接的界线勾画出来。再则是拼接影像图,当选好拼接线后,由专用软件沿着拼接线的轨迹自动进行拼接处理。最后是拼接后的检查,着重检查沿拼接线的接缝处是否存在着错位,若存在错位,还需要对拼接后的影像作进一步的修补。 

        在进行图像拼接时,必须注意以下三个问题: 一是拼接线要尽可能沿着道路、河流、田埂、空地、阴影等延伸,尽量将拼接线选择在两旁无高楼的区域。二是注意两幅相邻影像图在拼接处的高楼单中心投影倾向,要尽可能使拼接线两侧的楼房保持相似的倾向,同时也要防止在拼接后将某一侧楼房切掉一部分的情况。三是拼接线要尽可能避免穿越高架、桥梁、铁路等地物,假如必须要穿越高架、桥梁、铁路时,应尽量从衔接较好的地方或阴影区内穿过。 

        在进行错位修补时,必须遵循以下四个原则: 一要遵循客观性的原则,即在尊重原始影像的基础上,经过对影像错位的修补,使修补后的影像能客观地体现地物的原有面貌。二要遵循准确性的原则,即只对几何校正不准的影像部分作错位修补处理,而对几何校正准确的影像部分保持其原状不动。三要遵循整体性的原则,即无论是大尺度地物还是小尺度地物,只要拼接时在它的拼接处呈现错位,就要对整个地物作整体性地修补。四要遵循连锁性的原则,即对原有的影像错位作了锁定修补后,不要在其它地方再产生新的错位。

    ——————————————————————————————————————

    图像拼接技术主要包括两个关键环节即图像配准和图像融合对于图像融合部分,由于其耗时不太大,且现有的几种主要方法效果差别也不多,所以总体来说算法上比较成熟。而图像配准部分是整个图像拼接技术的核心部分,它直接关系到图像拼接算法的成功率和运行速度,因此配准算法的研究是多年来研究的重点。

           目前的图像配准算法基本上可以分为两类:基于频域的方法(相位相关方法)和基于时域的方法。
           
           相位相关法最早是由Kuglin和Hines在1975年提出的,并且证明在纯二维平移的情形下,拼接精度可以达到1个像素,多用于航空照片和卫星遥感图像的配准等领域。该方法对拼接的图像进行快速傅立叶变换,将两幅待配准图像变换到频域,然后通过它们的互功率谱直接计算出两幅图像间的平移矢量,从而实现图像的配准。由于其具有简单而精确的特点,后来成为最有前途的图像配准算法之一。但是相位相关方法一般需要比较大的重叠比例(通常要求配准图像之间有50%的重叠比例),如果重叠比例较小,则容易造成平移矢量的错误估计,从而较难实现图像的配准。

     

           基于时域的方法又可具体分为基于特征的方法和基于区域的方法。基于特征的方法首先找出两幅图像中的特征点(如边界点、拐点),并确定图像间特征点的对应关系,然后利用这种对应关系找到两幅图像间的变换关系。这一类方法不直接利用图像的灰度信息,因而对光线变化不敏感,但对特征点对应关系的精确程度依赖很大。这一类方法采用的思想较为直观,目前大部分的图像配准算法都可以归为这一类。基于区域的方法是以一幅图像重叠区域中的一块作为模板,在另一幅图像中搜索与此模板最相似的匹配块,这种算法精度较高,但计算量过大。

     

           按照匹配算法的具体实现又可以分为直接法和搜索法两大类,直接法主要包括变换优化法,它首先建立两幅待拼接图像间的变换模型,然后采用非线性迭代最小化算法直接计算出模型的变换参数,从而确定图像的配准位置。该算法效果较好,收敛速度较快,但是它要达到过程的收敛要求有较好的初始估计,如果初始估计不好,则会造成图像拼接的失败。搜索法主要是以一幅图像中的某些特征为依据,在另一幅图像中搜索最佳配准位置,常用的有比值匹配法,块匹配法和网格匹配法。比值匹配法是从一幅图像的重叠区域中部分相邻的两列上取出部分像素,然后以它们的比值作模板,在另一幅图像中搜索最佳匹配。这种算法计算量较小,但精度较低;块匹配法则是以一幅图像重叠区域中的一块作为模板,在另一幅图像中搜索与此模板最相似的匹配块,这种算法精度较高,但计算量过大;网格匹配法减小了块匹配法的计算量,它首先要进行粗匹配,每次水平或垂直移动一个步长,记录最佳匹配位置,然后在此位置附近进行精确匹配,每次步长减半,然后循环此过程直至步长减为0。这种算法较前两种运算量都有所减小,但在实际应用中仍然偏大,而且粗匹配时如果步长取的太大,很可能会造成较大的粗匹配误差,从而很难实现精确匹配。

    ——————————————————————————————————————

    网上搜的都是一行代码Stitcher::Status status = stitcher.stitch(imgs, pano);就出来的傻瓜拼接,连opencv基本的包都没用。

    自己好歹用了下基本的包实现了下。

    鲁棒性不太好,图片少的时候没事,图片一多就出现了内存错误和木有特征点的错误。

    #include <iostream>
    #include <fstream>
    #include <io.h>
    #include <stdio.h>
    #include <stdlib.h>
    #include <time.h>
    #include "opencv2/highgui/highgui.hpp"
    #include "opencv2/stitching/stitcher.hpp"
    #include "opencv2/stitching/warpers.hpp"

    #define MaxPics 30
    //F:\Program Files\opencv\modules\stitching\include\opencv2\stitching
    using namespace std;
    using namespace cv;
    bool try_use_gpu = false;
    int filenum;
    vector<Mat> imgs;
    string result_name = "result.jpg";
    //void printUsage();
    //int parseCmdArgs(int argc, char** argv);
    int readDir(char *path,char fnames[MaxPics][100]);

    int readDir(char *path,char fnames[MaxPics][100]) 

    struct _finddata_t c_file;
    intptr_t hFile;
    filenum=0;

    char fullpath[100];

    sprintf(fullpath,"%s\\*.jpg",path);
    // Find first .c file in current directory 
    if( (hFile = _findfirst( fullpath, &c_file )) == -1L )
    printf( "No *.jpg files in current directory %s!\n",fullpath);
    else
    {
    printf( "Listing of .jpg files in %s directory\n\n",path);
    printf( "RDO HID SYS ARC FILE SIZE\n", ' ' );
    printf( "--- --- --- --- ---- ----\n", ' ' );
    do {
    printf( ( c_file.attrib & _A_RDONLY ) ? " Y " : " N " );
    printf( ( c_file.attrib & _A_SYSTEM ) ? " Y " : " N " );
    printf( ( c_file.attrib & _A_HIDDEN ) ? " Y " : " N " );
    printf( ( c_file.attrib & _A_ARCH ) ? " Y " : " N " );
    //ctime_s( buffer, _countof(buffer), &c_file.time_write );
    printf( " %-15s %9ld\n",c_file.name, c_file.size );
    sprintf(fnames[filenum],"%s\\%s",path,c_file.name);
    filenum++;
    } while( _findnext( hFile, &c_file ) == 0 );
    _findclose( hFile );
    }
    return filenum;

    int main(int argc, char* argv[])
    {
    clock_t start,finish;
    double totaltime;
    start=clock();


    argv[1]="pic";
    char fdir[MaxPics][100];
    filenum = readDir(argv[1],fdir);
    Mat img, pano;
    for(int i=0 ; i<filenum ;i++){
    img = imread(fdir[i]);
    imgs.push_back(img);
    }
    Stitcher stitcher = Stitcher::createDefault(try_use_gpu);
    stitcher.setRegistrationResol(0.1);//为了加速,我选0.1,默认是0.6,最大值1最慢,此方法用于特征点检测阶段,如果找不到特征点,调高吧
    //stitcher.setSeamEstimationResol(0.1);//默认是0.1
    //stitcher.setCompositingResol(-1);//默认是-1,用于特征点检测阶段,找不到特征点的话,改-1
    stitcher.setPanoConfidenceThresh(1);//默认是1,见过有设0.6和0.4的
    stitcher.setWaveCorrection(false);//默认是true,为加速选false,表示跳过WaveCorrection步骤
    //stitcher.setWaveCorrectKind(detail::WAVE_CORRECT_HORIZ);//还可以选detail::WAVE_CORRECT_VERT ,波段修正(wave correction)功能(水平方向/垂直方向修正)。因为setWaveCorrection设的false,此语句没用

    //找特征点surf算法,此算法计算量大,但对刚体运动、缩放、环境影响等情况下较为稳定
    detail::SurfFeaturesFinder *featureFinder = new detail::SurfFeaturesFinder();
    stitcher.setFeaturesFinder(featureFinder);

    //找特征点ORB算法,但是发现草地这组图,这个算法不能完成拼接
    //detail::OrbFeaturesFinder *featureFinder = new detail::OrbFeaturesFinder();
    //stitcher.setFeaturesFinder(featureFinder);

    //Features matcher which finds two best matches for each feature and leaves the best one only if the ratio between descriptor distances is greater than the threshold match_conf.
    detail::BestOf2NearestMatcher *matcher = new detail::BestOf2NearestMatcher(false, 0.5f/*=match_conf默认是0.65,我选0.8,选太大了就没特征点啦,0.8都失败了*/);
    stitcher.setFeaturesMatcher(matcher);

    // Rotation Estimation,It takes features of all images, pairwise matches between all images and estimates rotations of all cameras.
    //Implementation of the camera parameters refinement algorithm which minimizes sum of the distances between the rays passing through the camera center and a feature,这个耗时短
    stitcher.setBundleAdjuster(new detail::BundleAdjusterRay());
    //Implementation of the camera parameters refinement algorithm which minimizes sum of the reprojection error squares.
    //stitcher.setBundleAdjuster(new detail::BundleAdjusterReproj());

    //Seam Estimation
    //Minimum graph cut-based seam estimator
    //stitcher.setSeamFinder(new detail::GraphCutSeamFinder(detail::GraphCutSeamFinderBase::COST_COLOR));//默认就是这个
    //stitcher.setSeamFinder(new detail::GraphCutSeamFinder(detail::GraphCutSeamFinderBase::COST_COLOR_GRAD));//GraphCutSeamFinder的第二种形式
    //啥SeamFinder也不用,Stub seam estimator which does nothing.
    stitcher.setSeamFinder(new detail::NoSeamFinder);
    //Voronoi diagram-based seam estimator.
    //stitcher.setSeamFinder(new detail::VoronoiSeamFinder);

    //exposure compensators曝光补偿
    //stitcher.setExposureCompensator(new detail::BlocksGainCompensator());//默认的就是这个
    //不要曝光补偿
    stitcher.setExposureCompensator(new detail::NoExposureCompensator());
    //Exposure compensator which tries to remove exposure related artifacts by adjusting image intensities
    //stitcher.setExposureCompensator(new detail::detail::GainCompensator());
    //Exposure compensator which tries to remove exposure related artifacts by adjusting image block intensities 
    //stitcher.setExposureCompensator(new detail::detail::BlocksGainCompensator()); 

    //Image Blenders
    //Blender which uses multi-band blending algorithm 
    //stitcher.setBlender(new detail::MultiBandBlender(try_use_gpu));//默认的是这个
    //Simple blender which mixes images at its borders
    stitcher.setBlender(new detail::FeatherBlender());//这个简单,耗时少

    //柱面?球面OR平面?默认为球面
    PlaneWarper* cw = new PlaneWarper();
    //SphericalWarper* cw = new SphericalWarper();
    //CylindricalWarper* cw = new CylindricalWarper();
    stitcher.setWarper(cw);

    Stitcher::Status status = stitcher.estimateTransform(imgs);
    if (status != Stitcher::OK)
    {
    cout << "Can't stitch images, error code = " << int(status) << endl;
    return -1;
    }
    status = stitcher.composePanorama(pano);
    if (status != Stitcher::OK)
    {
    cout << "Can't stitch images, error code = " << int(status) << endl;
    return -1;
    }
    cout<<"程序开始";
    imwrite(result_name, pano);
    finish=clock();
    totaltime=(double)(finish-start)/CLOCKS_PER_SEC;
    cout<<"\n此程序的运行时间为"<<totaltime<<"秒!"<<endl;
    return 0;
    }

    至于其他的拼接方式,http://academy.nearsoft.com/project-updates/makingapanoramapicture写的挺好。我引用之,并总结如下

    1,Stitcher class from OpenCV library

    • Pros
      • Extremely simple. Build panorama in just one command.
      • It performs from finding features through blending in the selected projection by itself.
    • Cons
      • Uses SURF algorithm (patented, not allowed for commercial use).
      • Slow at searching for features, running on a smartphone

    2,BoofCV

    • Pros
      • Still easy to use.
      • Java native.
    • Cons
      • It blends without warping into the spherical projection.
      • Still uses SURF.

    当然自己实现最好,我的上述代码就是第一种了

    ——————————————————————————————————————
    opencv2.4.0以上的版本提供了stitcher类,可以很方便的实现几幅图像的拼接,关于这个类详细的介绍,可以参考文档: http://docs.opencv.org/2.4.2/modules/stitching/doc/high_level.html?highlight=stitcher#stitcher
        该类主要用的成员函数有createDefault,用于创建缺省参数的stitcher;estimatedTransform,用于匹配图像、分析摄像头旋转角度;composePanorama,生成最后的拼接图像。文档中提示如果对stitching的整过过程不熟悉的话,最好不要使用以上两个函数,直接使用stitch就行了。关于图像拼接的过程,涉及到特征点的提取、特征点匹配、图像融合等等比较复杂的过程,可以参考相关论文和期刊。
        在安装文件下,提供了图像拼接的例子:C:\opencv2.4.2\opencv\samples\cpp\stitching.cpp  
    配置好后直接运行就可以了:

    #include "stdafx.h"
    #include <iostream>
    #include <fstream>
    #include "opencv2/highgui/highgui.hpp"
    #include "opencv2/stitching/stitcher.hpp"


    using namespace std;
    using namespace cv;

    bool try_use_gpu = false;
    vector<Mat> imgs;
    string result_name = "result.jpg";

    //void printUsage();
    //int parseCmdArgs(int argc, char** argv);

    int main(int argc, char* argv[])
    {

    Mat img=imread("1.jpg");
    imgs.push_back(img);
    img=imread("2.jpg");
    imgs.push_back(img);
    img=imread("3.jpg");
    imgs.push_back(img);

    Mat pano;
    Stitcher stitcher = Stitcher::createDefault(try_use_gpu);
    Stitcher::Status status = stitcher.stitch(imgs, pano);

    if (status != Stitcher::OK)
    {
    cout << "Can't stitch images, error code = " << int(status) << endl;
    return -1;
    }

    imwrite(result_name, pano);
    return 0;
    }

    最终结果:
    示例程序047--用opencv的stitcher类实现图像拼接

    示例程序047--用opencv的stitcher类实现图像拼接

    示例程序047--用opencv的stitcher类实现图像拼接

    示例程序047--用opencv的stitcher类实现图像拼接
    ——————————————————————————————————————
    第一章  绪论

      1.1 图像拼接技术的研究背景及研究意义

      图像拼接(image mosaic)是一个日益流行的研究领域,他已经成为照相绘图学、计算机视觉、图像处理和计算机图形学研究中的热点。图像拼接解决的问题一般式,通过对齐一系列空间重叠的图像,构成一个无缝的、高清晰的图像,它具有比单个图像更高的分辨率和更大的视野。

      早期的图像拼接研究一直用于照相绘图学,主要是对大量航拍或卫星的图像的整合。近年来随着图像拼接技术的研究和发展,它使基于图像的绘制(IBR)成为结合两个互补领域——计算机视觉和计算机图形学的坚决焦点,在计算机视觉领域中,图像拼接成为对可视化场景描述(Visual Scene Representaions)的主要研究方法:在计算机形学中,现实世界的图像过去一直用于环境贴图,即合成静态的背景和增加合成物体真实感的贴图,图像拼接可以使IBR从一系列真是图像中快速绘制具有真实感的新视图。

      在军事领域网的夜视成像技术中,无论夜视微光还是红外成像设备都会由于摄像器材的限制而无法拍摄视野宽阔的图片,更不用说360 度的环形图片了。但是在实际应用中,很多时候需要将360 度所拍摄的很多张图片合成一张图片,从而可以使观察者可以观察到周围的全部情况。使用图像拼接技术,在根据拍摄设备和周围景物的情况进行分析后,就可以将通过转动的拍摄器材拍摄的涵盖周围360 度景物的多幅图像进行拼接,从而实时地得到超大视角甚至是360 度角的全景图像。这在红外预警中起到了很大的作用。

      微小型履带式移动机器人项目中,单目视觉不能满足机器人的视觉导航需要,并且单目视觉机器人的视野范围明显小于双目视觉机器人的视野。利用图像拼接技术,拼接机器人双目采集的图像,可以增大机器人的视野,给机器人的视觉导航提供方便。在虚拟现实领域中,人们可以利用图像拼接技术来得到宽视角的图像或360 度全景图像,用来虚拟实际场景。这种基于全景图的虚拟现实系统,通过全景图的深度信息抽取,恢复场景的三维信息,进而建立三维模型。这个系统允许用户在虚拟环境中的一点作水平环视以及一定范围内的俯视和仰视,同时允许在环视的过程中动态地改变焦距。这样的全景图像相当于人站在原地环顾四周时看到的情形。在医学图像处理方面,显微镜或超声波的视野较小,医师无法通过一幅图像进行诊视,同时对于大目标图像的数据测量也需要把不完整的图像拼接为一个整体。所以把相邻的各幅图像拼接起来是实现远程数据测量和远程会诊的关键环节圆。在遥感技术领域中,利用图像拼接技术中的图像配准技术可以对来自同一区域的两幅或多幅图像进行比较,也可以利用图像拼接技术将遥感卫星拍摄到的有失真地面图像拼接成比较准确的完整图像,作为进一步研究的依据。

    从以上方面可以看出,图像拼接技术的应用前景十分广阔,深入研究图像拼接技术有着很重要的意义

      1.2图像拼接算法的分类

      图像拼接作为这些年来图像研究方面的重点之一,国内外研究人员也提出了很多拼接算法。图像拼接的质量,主要依赖图像的配准程度,因此图像的配准是拼接算法的核心和关键。根据图像匹配方法的不同仁阔,一般可以将图像拼接算法分为以下两个类型:
      (1) 基于区域相关的拼接算法。
          这是最为传统和最普遍的算法。基于区域的配准方法是从待拼接图像的灰度值出发,对待配准图像中一块区域与参考图像中的相同尺寸的区域使用最小二乘法或者其它数学方法计算其灰度值的差异,对此差异比较后来判断待拼接图像重叠区域的相似程度,由此得到待拼接图像重叠区域的范围和位置,从而实现图像拼接。也可以通过FFT 变换将图像由时域变换到频域,然后再进行配准。对位移量比较大的图像,可以先校正图像的旋转,然后建立两幅图像之间的映射关系。
    当以两块区域像素点灰度值的差别作为判别标准时,最简单的一种方法是直接把各点灰度的差值累计起来。这种办法效果不是很好,常常由于亮度、对比度的变化及其它原因导致拼接失败。另一种方法是计算两块区域的对应像素点灰度值的相关系数,相关系数越大,则两块图像的匹配程度越高。该方法的拼接效果要好一些,成功率有所提高。
     (2) 基于特征相关的拼接算法。
         基于特征的配准方法不是直接利用图像的像素值,而是通过像素导出图像的特征,然后以图像特征为标准,对图像重叠部分的对应特征区域进行搜索匹配,该类拼接算法有比较高的健壮性和鲁棒性。
    基于特征的配准方法有两个过程:特征抽取和特征配准。首先从两幅图像中提取灰度变化明显的点、线、区域等特征形成特征集冈。然后在两幅图像对应的特征集中利用特征匹配算法尽可能地将存在对应关系的特征对选择出来。一系列的图像分割技术都被用到特征的抽取和边界检测上。如canny 算子、拉普拉斯高斯算子、区域生长。抽取出来的空间特征有闭合的边界、开边界、交叉线以及其他特征。特征匹配的算法有:交叉相关、距离变换、动态编程、结构匹配、链码相关等算法。

      1.3本文的主要工作和组织结构

    本文的主要工作:
    (1) 总结了前人在图像拼接方面的技术发展历程和研究成果。

    (2) 学习和研究了前人的图像配准算法。

    (3) 学习和研究了常用的图像融合算法。

    (4) 用matlab实现本文中的图像拼接算法

    (5) 总结了图像拼接中还存在的问题,对图像拼接的发展方向和应用前景进行展望。

      本文的组织结构:

      第一章主要对图像拼接技术作了整体的概述,介绍了图像拼接的研究背景和应用前景,以及图像拼接技术的大致过程、图像拼接算法的分类和其技术难点。第二章主要介绍讨论了图像预处理中的两个步骤,即图像的几何校正和噪声点的抑制。第三章主要介绍讨论了图像配准的多种算法。第四章主要介绍讨论了图像融合的一些算法。第五章主要介绍图像拼接软件实现本文的算法。第六章主要对图像拼接中还存在的问题进行总结,以及对图像拼接的发展进行展望。

      1.4 本章小结

      本章主要对图像拼接技术作了整体的概述,介绍了图像拼接的研究背景和应用前景,以图像拼接算法的分类和其技术难点,并且对全文研究内容进行了总体介绍。
      第二章  图像拼接的基础理论及图像预处理

      2.1图像拼接

        图像拼接技术主要有三个主要步骤:图像预处理、图像配准、图像融合与边界平滑,

    如图。

    图像拼接算法及实现

        图像拼接技术主要分为三个主要步骤:图像预处理、图像配准、图像融合与边界平滑,图像预处理主要指对图像进行几何畸变校正和噪声点的抑制等,让参考图像和待拼接图像不存在明显的几何畸变。在图像质量不理想的情况下进行图像拼接,如果不经过图像预处理,很容易造成一些误匹配。图像预处理主要是为下一步图像配准做准备,让图像质量能够满足图像配准的要求。图像配准主要指对参考图像和待拼接图像中的匹配信息进行提取,在提取出的信息中寻找最佳的匹配,完成图像间的对齐。图像拼接的成功与否主要是图像的配准。待拼接的图像之间,可能存在平移、旋转、缩放等多种变换或者大面积的同色区域等很难匹配的情况,一个好的图像配准算法应该能够在各种情况下准确找到图像间的对应信息,将图像对齐。图像融合指在完成图像匹配以后,对图像进行缝合,并对缝合的边界进行平滑处理,让缝合自然过渡。由于任何两幅相邻图像在采集条件上都不可能做到完全相同,因此,对于一些本应该相同的图像特性,如图像的光照特性等,在两幅图像中就不会表现的完全一样。图像拼接缝隙就是从一幅图像的图像区域过渡到另一幅图像的图像区域时,由于图像中的某些相关特性发生了跃变而产生的。图像融合就是为了让图像间的拼接缝隙不明显,拼接更自然

      2.2 图像的获取方式

      图像拼接技术原理是根据图像重叠部分将多张衔接的图像拼合成一张高分辨率全景图 。这些有重叠部分的图像一般由两种方法获得 : 一种是固定照相机的转轴 ,然后绕轴旋转所拍摄的照片 ;另一种是固定照相机的光心 ,水平摇动镜头所拍摄的照片。其中 ,前者主要用于远景或遥感图像的获取 ,后者主要用于显微图像的获取 ,它们共同的特点就是获得有重叠的二维图像。

      2.3 图像的预处理

      2.3.1 图像的校正

      当照相系统的镜头或者照相装置没有正对着待拍摄的景物时候,那么拍摄到的景物图像就会产生一定的变形。这是几何畸变最常见的情况。另外,由于光学成像系统或电子扫描系统的限制而产生的枕形或桶形失真,也是几何畸变的典型情况。几何畸变会给图像拼接造成很大的问题,原本在两幅图像中相同的物体会因为畸变而变得不匹配,这会给图像的配准带来很大的问题。因此,解决几何畸变的问题显得很重要。

      图象校正的基本思路是,根据图像失真原因,建立相应的数学模型,从被污染或畸变的图象信号中提取所需要的信息,沿着使图象失真的逆过程恢复图象本来面貌。实际的复原过程是设计一个滤波器,使其能从失真图象中计算得到真实图象的估值,使其根据预先规定的误差准则,最大程度地接近真实图象。

      2.3.2 图像噪声的抑制

      图像噪声可以理解为妨碍人的视觉感知,或妨碍系统传感器对所接受图像源信息进行理解或分析的各种因素,也可以理解成真实信号与理想信号之间存在的偏差。一般来说,噪声是不可预测的随机信号,通常采用概率统计的方法对其进行分析。噪声对图像处理十分重要,它影响图像处理的各个环节,特别在图像的输入、采集中的噪声抑制是十分关键的问题。若输入伴有较大的噪声,必然影响图像拼接的全过程及输出的结果。根据噪声的来源,大致可以分为外部噪声和内部噪声;从统计数学的观点来定义噪声,可以分为平稳噪声和非平稳噪声。各种类型的噪声反映在图像画面上,大致可以分为两种类型。一是噪声的幅值基本相同,但是噪声出现的位置是随机的,一般称这类噪声为椒盐噪声。另一种是每一点都存在噪声,但噪声的幅值是随机分布的,从噪声幅值大小的分布统计来看,其密度函数有高斯型、瑞利型,分别成为高斯噪声和瑞利噪声,又如频谱均匀分布的噪声称为白噪声等。

      1.均值滤波

        所谓均值滤波实际上就是用均值替代原图像中的各个像素值。均值滤波的方法是,对将处理的当前像素,选择一个模板,该模板为其邻近的若干像素组成,用模板中像素的均值来替代原像素的值。如图2.4所示,序号为0是当前像素,序号为1至8是邻近像素。求模板中所有像素的均值,再把该均值赋予当前像素点((x, y),作为处理后图像在该点上的灰度g(x,y),即

     g(x,y)=图像拼接算法及实现 图像拼接算法及实现                                                      (2-2-2-1)

    其中,s为模板,M为该模板中包含像素的总个数。

    图像拼接算法及实现

    图2.2.2.1模板示意图

      2.中值滤波

      中值滤波是基于排序统计理论的一种能有效抑制噪声的非线性信号处理技术。它的核心算法是将模板中的数据进行排序,这样,如果一个亮点(暗点)的噪声,就会在排序过程中被排在数据序列的最右侧或者最左侧,因此,最终选择的数据序列中见位置上的值一般不是噪声点值,由此便可以达到抑制噪声的目的。
    取某种结构的二维滑动模板,将模板内像素按照像素值的大小进行排序,生成单调上升(或下降)的二维数据序列。二维德中值滤波输出为
                                     ( 2-2-2-2 )                            

    其中,f(x,y),g (x,y)分别为原图像和处理后的图像,w二维模板,k ,l为模板的长宽,Med 为取中间值操作,模板通常为3  3 、5 5 区域,也可以有不同形状,如线状、圆形、十字形、圆环形。

      2.4 本章小结

      本章主要介绍了图像几何畸变校正和图像噪声抑制两种图像预处理.

      第三章  图像配准算法

      3.1 图像配准的概念

        图像配准简而言之就是图像之间的对齐。图像配准定义为:对从不同传感器或不同时间或不同角度所获得的两幅或多幅图像进行最佳匹配的处理过程。为了更清楚图像配准的任务,我们将图像配准问题用更精确的数学语言描述出来。配准可以用描述为如下的问题:

    给定同一景物的从不同的视角或在不同的时间获取的两个图像I ,I 和两个图像间的相似度量S(I ,I ),找出I ,I 中的同名点,确定图像间的最优变换T,使得S(T(I ),I )达到最大值。

        图像配准总是相对于多幅图像来讲的,在实际工作中,通常取其中的一幅图像作为配准的基准,称它为参考图,另一幅图像,为搜索图。图像配准的一般做法是,首先在参考图上选取以某一目标点为中心的图像子块,并称它为图像配准的模板,然后让模板在搜索图上有秩序地移动,每移到一个位置,把模板与搜索图中的对应部分进行相关比较,直到找到配准位置为止。

        如果在模板的范围内,同一目标的两幅图像完全相同,那么完成图像配准并不困难。然而,实际上图像配准中所遇到的同一目标的两幅图像常常是在不同条件下获得的,如不同的成像时间、不同的成像位置、甚至不同的成像系统等,再加上成像中各种噪声的影响,使同一目标的两幅图像不可能完全相同,只能做到某种程度的相似,因此图像配准是一个相当复杂的技术过程。
      3.2 基于区域的配准

      3.2.1 逐一比较法

    设搜索图为s待配准模板为T,如图3.1所示,S大小为M N,T大小为U V,如图所示。

    图像拼接算法及实现

                  图3.1搜索图S与模板T示意图

     逐一比较法的配准思想是:

        在搜索图S中以某点为基点(i,j),截取一个与模板T大小一样的分块图像,这样的基点有(M-U+1) 图像拼接算法及实现(N-V+1)个,配准的目标就是在(M-U+1) 图像拼接算法及实现(N-V+1)个分块图像中找一个与待配准图像最相似的图像,这样得到的基准点就是最佳配准点。

        设模板T在搜索图s上移动,模板覆盖下的那块搜索图叫子图S ,(i,j)为这块子图的左上角点在S图中的坐标,叫做参考点。然后比较T和S 图像拼接算法及实现的内容。若两者一致,则T和S图像拼接算法及实现 之差为零。在现实图像中,两幅图像完全一致是很少见的,一般的判断是在满足一定条件下,T和S图像拼接算法及实现 之差最小。

    根据以上原理,可采用下列两种测度之一来衡量T和S 的相似程度。D(i,j)的值越小,则该窗口越匹配。

              D(i,j)= 图像拼接算法及实现 图像拼接算法及实现  [S图像拼接算法及实现 (m,n)-T(m,n)]图像拼接算法及实现              (3-1)

              D(i,j)=图像拼接算法及实现 图像拼接算法及实现  图像拼接算法及实现 [S图像拼接算法及实现 (m,n)-T(m,n)图像拼接算法及实现              (3-2)

    或者利用归一化相关函数。将式(3-1)展开可得:

    D(i,j)=图像拼接算法及实现 图像拼接算法及实现  [S (m,n)] -2 图像拼接算法及实现 图像拼接算法及实现 S (m,n)*T(m,n)+图像拼接算法及实现 图像拼接算法及实现  [T(m,n)]  (3-3)

        式中等号右边第三项表示模板总能量,是一常数,与(i,j)无关;第一项是与模板匹配区域的能量,它随((i,j)的改变而改变,当T和S 匹配时的取最大值。因此相

    关函数为:

    R(i,j)=图像拼接算法及实现                                             (3-4)

    当R(i,j)越大时,D(i,j)越小,归一化后为:

     R(i,j)=图像拼接算法及实现                      (3-5)

    根据Cauchy-Schwarz不等式可知式(3-5)中0图像拼接算法及实现 R(i,j)图像拼接算法及实现 1,并且仅当值S (m, n)/T (m, n)=常数时,R(i,j)取极大值。

        该算法的优点:

        (1)算法思路比较简单,容易理解,易于编程实现。

        (2)选用的模板越大,包含的信息就越多,匹配结果的可信度也会提高,同时能够对参考图像进行全面的扫描。

        该算法的缺点:

        (1)很难选择待配准图像分块。因为一个如果分块选择的不正确,缺少信息量,则不容易正确的匹配,即发生伪匹配。同时,如果分块过大则降低匹配速度,如果分块过小则容易降低匹配精度。

        (2)对图像的旋转变形不能很好的处理。算法本身只是把待配准图像分块在标准参考图像中移动比较,选择一个最相似的匹配块,但是并不能够对图像的旋转变形进行处理,因此对照片的拍摄有严格的要求。

      3.2.2 分层比较法

      图像处理的塔形(或称金字塔:Pyramid)分解方法是由Burt和Adelson首先提出的,其早期主要用于图像的压缩处理及机器人的视觉特性研究。该方法把原始图像分解成许多不同空间分辨率的子图像,高分辨率(尺寸较大)的子图像放在下层,低分辨率(尺寸较小)的图像放在上层,从而形成一个金字塔形状。

        在逐一比较法的思想上,为减少运算量,引入了塔形处理的思想,提出了分层比较法。利用图像的塔形分解,可以分析图像中不同大小的物体。同时,通过对低分辨率、尺寸较小的上层进行分析所得到的信息还可以用来指导对高分辨率、尺寸较大的下层进行分析,从而大大简化分析和计算。在搜索过程中,首先进行粗略匹配,每次水平或垂直移动一个步长,计算对应像素点灰度差的平方和,记录最小值的网格位置。其次,以此位置为中心进行精确匹配。每次步长减半,搜索当前最小值,循环这个过程,直到步长为零,最后确定出最佳匹配位置。

        算法的具体实现步骤如下:

        (1)将待匹配的两幅图像中2 图像拼接算法及实现 2邻域内的像素点的像素值分别取平均,作为这一区域(2 图像拼接算法及实现 2)像素值,得到分辨率低一级的图像。然后,将此分辨率低一级的图像再作同样的处理,也就是将低一级的图像4 图像拼接算法及实现 4邻域内的像素点的像素值分别取平均,作为这一区域(4图像拼接算法及实现  4)点的像素值,得到分辨率更低一级的图像。依次处理,得到一组分辨率依次降低的图像。

        (2)从待匹配的两幅图像中分辨率最低的开始进行匹配搜索,由于这两幅图像像素点的数目少,图像信息也被消除一部分,因此,此匹配位置是不精确的。所以,在分辨率更高一级的图像中搜索时,应该在上一次匹配位置的附近进行搜索。依次进行下去,直到在原始图像中寻找到精确的匹配位置。

        算法的优点:

          (1)该算法思路简单,容易理解,易于编程实现。

          (2)该算法的搜索空间比逐一比较要少,在运算速度较逐一比较法有所提高。

        算法的缺点:

          (1)算法的精度不高。在是在粗略匹配过程中,移动的步长较大,很有可能将第一幅图像上所取的网格划分开,这样将造成匹配中无法取出与第一幅图像网格完全匹配的最佳网格,很难达到精确匹配。

          (2)对图像的旋转变形仍然不能很好的处理。与逐一比较法一样,该算法只是对其运算速度有所改进,让搜索空间变小,并无本质变化,因此对图像的旋转变形并不能进行相应处理。

      3.2.3 相位相关法

      相位相关度法是基于频域的配准常用算法。它将图像由空域变换到频域以后再进行配准。该算法利用了互功率谱中的相位信息进行图像配准,对图像间的亮度变化不敏感,具有一定的抗干扰能力,而且所获得的相关峰尖锐突出,位移检测范围大,具有较高的匹配精度。

        相位相关度法思想是利用傅立叶变换的位移性质,对于两幅数字图像s,t,其对应的傅立叶变换为S,T,即:

    S=F{s}=图像拼接算法及实现 e图像拼接算法及实现    T=F{t}=图像拼接算法及实现 e图像拼接算法及实现              (3-6)

    若图像s,t相差一个平移量(x ,y ),即有:

          s(x,y) = t(x-x ,y-y )                               (3-7)

        根据傅立叶变换的位移性质,上式的傅立叶变换为:

         S(图像拼接算法及实现 )=e图像拼接算法及实现 T(图像拼接算法及实现 )                        (3-8)

        也就是说,这两幅图像在频域中具有相同的幅值,只是相位不同,他们之间的相位差可以等效的表示为互功率谱的相位。两幅图的互功率谱为:

     图像拼接算法及实现=e 图像拼接算法及实现                       (3-9)

        其中*为共扼符号, 表示频谱幅度。通过对互功率谱式(3-9)进行傅立叶逆变换,在((x,y)空间的(x ,y ),即位移处,将形成一个单位脉冲函数 ,脉冲位置即为两幅被配准图像间的相对平移量x 和y

        式(3-9)表明,互功率谱的相位等价于图像间的相位差,故该方法称作相位相关法。

        相位相关度法的优点:

          (1)该算法简单速度快,因此经常被采用。对于其核心技术傅立叶变换,现在己经出现了很多有关的快速算法,这使得该算法的快速性成为众多算法中的一大优势。另外,傅立叶变换的硬件实现也比其它算法容易。

          (2)该算法抗干扰能力强,对于亮度变化不敏感。

      相位相关度法的缺点:

        (1)该算法要求图像有50%左右的重叠区域,在图像重叠区域很小的时,算法的结果很难保证,容易造成误匹配。

        (2)由于Fourier变换依赖于自身的不变属性,所以该算法只适用于具有旋转、平移、比例缩放等变换的图像配准问题。对于任意变换模型,不能直接进行处理,而要使用控制点方法,控制点方法可以解决诸如多项式、局部变形等问题。

      3.3 基于特征的配准

      3.3.1 比值匹配法

      比值匹配法算法思路是利用图像中两列上的部分像素的比值作为模板,即在参考图像T的重叠区域中分别在两列上取出部分像素,用它们的比值作为模板,然后在搜索图S中搜索最佳的匹配。匹配的过程是在搜索图S中,由左至右依次从间距相同的两列上取出部分像素,并逐一计算其对应像素值比值;然后将这些比值依次与模板进行比较,其最小差值对应的列就是最佳匹配。这样在比较中只利用了一组数据,而这组数据利用了两列像素及其所包含的区域的信息。

        该算法的具体实现步骤如下:

        (1)在参考图像T中间隔为c个像素的距离上的两列像素中,各取m个像素,计算这m个像素的比值,将m个比值存入数组中,将其作为比较的模板。

        (2)从搜索图S中在同样相隔c个像素的距离上的两列,各取出m+n个像素,计算其比值,将m+n个比值存入数组。假定垂直错开距离不超过n个像素,多取的n个像素则可以解决图像垂直方向上的交错问题。

        (3)利用参考图像T中的比值模板在搜索图S中寻找相应的匹配。首先进行垂直方向上的比较,即记录下搜索图S中每个比值数组内的最佳匹配。再将每个数组的组内最佳匹配进行比较,即进行水平方向的比较,得到的最小值就认为是全局最佳匹配。此时全局最佳匹配即为图像间在水平方向上的偏移距离,该全局最佳匹配队应的组内最佳匹配即为图像间垂直方向上的偏移距离。

        比值匹配法的优点:

        (1)算法思路清晰简单,容易理解,实现起来比较方便。

        (2)在匹配计算的时候,计算量小,速度快。

        比值匹配法的缺点:

        (1)利用图像的特征信息太少。只利用了两条竖直的平行特征线段的像素的信息,没有能够充分利用了图像重叠区域的大部分特征信息。虽然算法提到,在搜索图S中由左至右依次从间距相同的两列上取出部分像素,计算其对应像素的比值,然后将这些比值依次与模版进行比较,好像是利用了搜索图S中的重叠区域的大部分图像信息,但在参考图像T中,只是任意选择了两条特征线,没有充分利用到参考图像T的重叠区域的特征信息。

        (2)对图片的采集提出了较高的要求。此算法对照片先进行垂直方向上的比较,然后再进行水平方向上的比较,这样可以解决上下较小的错开问题。在采集的时候只能使照相机在水平方向上移动。然而,有时候不可避免的照相机镜头会有小角度的旋转,使得拍摄出来的照片有一定的旋转,在这个算法中是无法解决的。而且对重叠区域无明显特征的图像,比较背景是海洋或者天空,这样在选取特征模版的时候存在很大的问题。由于照片中存在大块纹理相同的部分,所以与模版的差别就不大,这样有很多匹配点,很容易造成误匹配。

        (3)不易对两条特征线以及特征线之间的距离进行确定。算法中在参考图像T的重叠区域中取出两列像素上的部分像素,并没有给出选择的限制。然而在利用拼接算法实现自动拼接的时候,如果选取的特征线不是很恰当,那么这样的特征线算出来的模版就失去了作为模版的意义。同时,在确定特征线间距时,选的过大,则不能充分利用重叠区域的图像信息。选择的过小,则计算量太大。

     

    来源:http://www.ce86.com/a/computer/jsjyy123/200906/17-72995-2.html

    ——————————————————————————————————————

    基本介绍

    图像配准图像融合是图像拼接的两个关键技术。图像配准是图像融合的基础,而且图像配准算法的计算量一般非常大,因此图像拼接技术的发展很大程度上取决于图像配准技术的创新。早期的图像配准技术主要采用点匹配法,这类方法速度慢、精度低,而且常常需要人工选取初始匹配点,无法适应大数据量图像的融合。图像拼接的方法很多,不同的算法步骤会有一定差异,但大致的过程是相同的。一般来说,图像拼接主要包括以下五步:

    折叠图像预处理

    包括数字图像处理的基本操作(如去噪、边缘提取、直方图处理等)、建立图像的匹配模板以及对图像进行某种变换(如傅里叶变换、小波变换等)等操作。

    折叠图像配准

    就是采用一定的匹配策略,找出待拼接图像中的模板或特征点在参考图像中对应的位置,进而确定两幅图像之间的变换关系。

    折叠建立变换模型

    根据模板或者图像特征之间的对应关系,计算出数学模型中的各参数值,从而建立两幅图像的数学变换模型。

    折叠统一坐标变换

    根据建立的数学转换模型,将待拼接图像转换到参考图像的坐标系中,完成统一坐标变换。

    折叠融合重构

    将带拼接图像的重合区域进行融合得到拼接重构的平滑无缝全景图像。

    折叠编辑本段图像拼接的基本流程图。

    相邻图像的配准及拼接是全景图生成技术的关键,有关图像配准技术的研究至今已有很长的历史,其主要的方法有以下两种:基于两幅图像的亮度差最小的方法和基于特征的方法。其中使用较多的是基于特征模板匹配特征点的拼接方法。该方法允许待拼接的图像有一定的倾斜和变形,克服了获取图像时轴心必须一致的问题,同时允许相邻图像之间有一定色差。全景图的拼接主要包括以下4个步骤:图像的预拼接,即确定两幅相邻图像重合的较精确位置,为特征点的搜索奠定基础。特征点的提取,即在基本重合位置确定后,找到待匹配的特征点。图像矩阵变换及拼接,即根据匹配点建立图像的变换矩阵并实现图像的拼接。最后是图像的平滑处理。

    折叠编辑本段图像拼接技术分类

    图像拼接技术主要包括两个关键环节即图像配准和图像融合。对于图像融合部分,由于其耗时不太大,且现有的几种主要方法效果差别也不多,所以总体来说算法上比较成熟。而图像配准部分是整个图像拼接技术的核心部分,它直接关系到图像拼接算法的成功率和运行速度,因此配准算法的研究是多年来研究的重点。

    目前的图像配准算法基本上可以分为两类:基于频域的方法(相位相关方法)和基于时域的方法。

    折叠相位相关法

    相位相关法最早是由Kuglin和Hines在1975年提出的,并且证明在纯二维平移的情形下,拼接精度可以达到1个像素,多用于航空照片和卫星遥感图像的配准等领域。该方法对拼接的图像进行快速傅立叶变换,将两幅待配准图像变换到频域,然后通过它们的互功率谱直接计算出两幅图像间的平移矢量,从而实现图像的配准。由于其具有简单而精确的特点,后来成为最有前途的图像配准算法之一。但是相位相关方法一般需要比较大的重叠比例(通常要求配准图像之间有50%的重叠比例),如果重叠比例较小,则容易造成平移矢量的错误估计,从而较难实现图像的配准。

    折叠基于时域的方法

    基于时域的方法又可具体分为基于特征的方法和基于区域的方法。基于特征的方法首先找出两幅图像中的特征点(如边界点、拐点),并确定图像间特征点的对应关系,然后利用这种对应关系找到两幅图像间的变换关系。这一类方法不直接利用图像的灰度信息,因而对光线变化不敏感,但对特征点对应关系的精确程度依赖很大。这一类方法采用的思想较为直观,大部分的图像配准算法都可以归为这一类。基于区域的方法是以一幅图像重叠区域中的一块作为模板,在另一幅图像中搜索与此模板最相似的匹配块,这种算法精度较高,但计算量过大。

    按照匹配算法的具体实现又可以分为直接法和搜索法两大类,直接法主要包括变换优化法,它首先建立两幅待拼接图像间的变换模型,然后采用非线性迭代最小化算法直接计算出模型的变换参数,从而确定图像的配准位置。该算法效果较好,收敛速度较快,但是它要达到过程的收敛要求有较好的初始估计,如果初始估计不好,则会造成图像拼接的失败。搜索法主要是以一幅图像中的某些特征为依据,在另一幅图像中搜索最佳配准位置,常用的有比值匹配法,块匹配法和网格匹配法。比值匹配法是从一幅图像的重叠区域中部分相邻的两列上取出部分像素,然后以它们的比值作模板,在另一幅图像中搜索最佳匹配。这种算法计算量较小,但精度较低;块匹配法则是以一幅图像重叠区域中的一块作为模板,在另一幅图像中搜索与此模板最相似的匹配块,这种算法精度较高,但计算量过大;网格匹配法减小了块匹配法的计算量,它首先要进行粗匹配,每次水平或垂直移动一个步长,记录最佳匹配位置,然后在此位置附近进行精确匹配,每次步长减半,然后循环此过程直至步长减为0。这种算法较前两种运算量都有所减小,但在实际应用中仍然偏大,而且粗匹配时如果步长取的太大,很可能会造成较大的粗匹配误差,从而很难实现精确匹配。

    ——————————————————————————————————————

    stitching是OpenCV2.4.0一个新模块,功能是实现图像拼接,所有的相关函数都被封装在Stitcher类当中。这个类当中我们可能用到的成员函数有createDefault、estimateTransform、composePanorama、stitch。其内部实现的过程是非常繁琐的,需要很多算法的支持,包括图像特征的寻找和匹配,摄像机的校准,图像的变形,曝光补偿和图像融合。但这些模块的接口、调用,强大的OpenCV都为我们搞定了,我们使用OpenCV做图像拼接,只需要调用createDefault函数生成默认的参数,再使用stitch函数进行拼接就ok了。就这么简单!estimateTransform和composePanorama函数都是比较高级的应用,如果各位对stitching的流程不是很清楚的话,还是慎用。

    实例也非常简单,下载链接哦:http://download.csdn.net/detail/yang_xian521/4321158

    输入原图(为了显示,我都压缩过):





    展开全文
  • 接opencv6.4-imgproc图像处理模块之直方图与模板 这部分的《opencv_tutorial》上都是直接上代码,没有原理部分的解释的。 十一、轮廓 1、图像中找轮廓 /// 转成灰度并模糊化降噪 cvtColor( src, src_...

    opencv6.4-imgproc图像处理模块之直方图与模板

    这部分的《opencv_tutorial》上都是直接上代码,没有原理部分的解释的。

    十一、轮廓

    1、图像中找轮廓


      /// 转成灰度并模糊化降噪
      cvtColor( src, src_gray, CV_BGR2GRAY );
      blur( src_gray, src_gray, Size(3,3) );
    Mat canny_output;//找到轮廓的图
      vector<vector<Point> > contours;//装载曲线的点
      vector<Vec4i> hierarchy;
    
      /// 用Canny算子检测边缘
      Canny( src_gray, canny_output, thresh, thresh*2, 3 );
      /// 寻找轮廓
      findContours( canny_output, contours, hierarchy, CV_RETR_TREE, CV_CHAIN_APPROX_SIMPLE, Point(0, 0) );
    
      /// 绘出轮廓
      Mat drawing = Mat::zeros( canny_output.size(), CV_8UC3 );
      for( int i = 0; i< contours.size(); i++ )//通过对contours.size(),就知道有几个分离的轮廓了
         {
           Scalar color = Scalar( rng.uniform(0, 255), rng.uniform(0,255), rng.uniform(0,255) );
           drawContours( drawing, contours, i, color, 2, 8, hierarchy, 0, Point() );//画出轮廓
         }
    
      /// 在窗体中显示结果
      namedWindow( "Contours", CV_WINDOW_AUTOSIZE );
      imshow( "Contours", drawing );
    查找轮廓的函数原型:void findContours(InputOutputArray image, OutputArrayOfArrays contours, OutputArray hierarchy, int mode, int method, Point offset=Point())
                                             void findContours(InputOutputArray image, OutputArrayOfArrays contours, int mode, int method, Point offset=Point())

    参数列表:image:一个8-bit的单通道原图像。非零的像素值被当作1来处理(类似bool类型)。0的像素值仍旧还是0,所以这个image是被当作二值处理的。你可以使用compare(),inRange(),threshold(),adaptiveThreshold(),Canny(),和其他的从一个灰度或者彩色图中生成一个二值图。该函数会在提出轮廓的时候修改image。如果mode等于CV_RETR_CCOMP或者CV_RETR_FLOODFILL,那么这个输入同样可以是标签的32-bit整数图(CV_32SC1)。

                      contours:检测到的轮廓,每个轮廓都是存储为点向量;

                      hierarchy:可选的输出向量,包含图像拓扑信息。它有着和轮廓数量一样多的元素。对于第i 个轮廓contours[i],元素hierarchy【i】【0】,hiearchy【i】【1】,hiearchy【i】【2】和hiearchy【i】【3】是在contours中基于0索引的,分别表示在同一个层次级别上的下一个轮廓上一个轮廓,第一个孩子轮廓父亲轮廓。如果轮廓 i 没有下一个、上一个、父亲或者嵌套的轮廓,对应的hierarchy【i】的元素就是负的。(这里其实就是个树结构,来进行索引不同的轮廓)

                      mode:轮廓索引模型

                     – CV_RETR_EXTERNAL     只检索最外面的轮廓。它会对所有的轮廓设置 hierarchy[i][2] = hierarchy[i][3] = -1 .
                     – CV_RETR_LIST                 不建立任何层次关系来进行索引所有的轮廓.
                     – CV_RETR_CCOMP       索引所有的轮廓并将它们组织成一个two-level的hierarchy(两个层次的层级关系):在顶层上,有着成分的外部边界;在第二层是孔洞的边界。如果有另一个轮廓在一个连接起来的成分内部,那么就将它放在顶层上。
                    – CV_RETR_TREE        检索所有的轮廓并重构一个全层次的嵌套轮廓结构,在contours.c的例子中你更可以看到全层次建立的代码。
                     method :轮廓逼近的方法

                       – CV_CHAIN_APPROX_NONE   存储所有轮廓点的绝对值。也就是说任何的轮廓的2子序列点 (x1,y1) 和 (x2,y2) 就是表示水平的,竖直的或者对角线邻居,也就是: max(abs(x1-x2),abs(y2-y1))==1;
                    – CV_CHAIN_APPROX_SIMPLE    压缩水平的,竖直的和对角线线段,只保留他们的终端。例如:对于一个up-right 的矩形轮廓是由4个点进行编码的。

                     – CV_CHAIN_APPROX_TC89_L1,  CV_CHAIN_APPROX_TC89_KCOS     使用 Teh-Chin 链逼近算法中主流之一的算法。

                    offset:可选的偏移量,通过这个值可以平移每一个轮廓。当轮廓是从一个图像的ROI中提取而你需要在整个图像上下文中分析的时候会变得很有用。

    该函数使用算法从二值图像中索引轮廓。轮廓对于形状分析和对象的检测识别是很有用的,在squares.c中有例子代码。

    Notes:image会被该函数所修改,同样的该函数不考虑图像的1像素边界(它会被0填充然后用来作为邻居分析),因此接触图像边界的轮廓会被修剪(clip,夹)。

    画出轮廓的函数原型:void drawContours(InputOutputArray image, InputArrayOfArrays contours, int contourIdx, const Scalar& color, int thickness=1, int lineType=LINE_8, InputArray hierarchy= noArray(), int maxLevel=INT_MAX, Point offset=Point() );

    参数列表:目标图像、输入的轮廓、轮廓的索引id、轮廓的颜色、粗度、线类型、可选的层次信息、最大级别、偏移量;

       第一个参数:目标图像,即画布;

       第二个参数:所有的输入轮廓,每个轮廓存储成点向量的形式;

       第三个参数:用来索引需要画的轮廓的参数,如果是负的,那么就画所有的轮廓;

       第四个参数:轮廓的颜色;

       第五个参数:画的轮廓的线的粗细程度,如果是负的(例如:thickness = CV_FILLED),轮廓内部也会画出来;

       第六个参数:线连接类型,line()可以有更详细的说明。

       第七个参数:可选的关于层次的信息。当你只想画某些轮廓的时候才是需要的(见maxLevel);

       第八个参数:画轮廓的最大等级。如果为0,那么只画指定的轮廓。如果为1,该函数会画轮廓及所有的嵌套轮廓。如果为2,该函数会画轮廓、所有的嵌套轮廓、所有的嵌套-to-嵌套的轮廓,以此类推。该参数当hierarchy可用的时候才被考虑;

       第九个参数:可选的轮廓平移参数。通过制定的平移量offset = (dx,dy)来平移所有的画的轮廓。

    该函数但thickness >=0的时候会画轮廓的外部线条,或者当thickness <0的时候会填充轮廓的框起来的区域。下面的代码是展示如何从二值图像中索引联通的成分,并且标记他们:

    Mat src;
         src = imread(name,0);
    Mat dst = Mat::zeros(src.rows, src.cols, CV_8UC3);
         src = src > 1;
    
    vector<vector<Point> > contours;
    vector<Vec4i> hierarchy;
    findContours( src, contours, hierarchy,RETR_CCOMP, CHAIN_APPROX_SIMPLE );
    //通过对top-level的轮廓进行迭代来画出每个联通的成分,使用随机颜色
    
    int idx = 0;
       for( ; idx >= 0; idx = hierarchy[idx][0] ){
         Scalar color( rand()&255, rand()&255, rand()&255 );
         drawContours( dst, contours, idx, color, FILLED, 8, hierarchy );
    }
         namedWindow( "Components", 1 );
         imshow( "Components", dst );
         waitKey(0);
    2、计算物体的凸包

    凸包就是将对象进行外部轮廓的包含,而且是凸图形的:


    /// 转成灰度图并进行模糊降噪
       cvtColor( src, src_gray, CV_BGR2GRAY );
       blur( src_gray, src_gray, Size(3,3) );
    Mat src_copy = src.clone();
       Mat threshold_output;
       vector<vector<Point> > contours;//存储轮廓的点集合
       vector<Vec4i> hierarchy;//构建轮廓的层次结构
    
       /// 对图像进行二值化
    int thresh = 100;
     int max_thresh = 255;
     RNG rng(12345);
    threshold( src_gray, threshold_output, thresh, 255, THRESH_BINARY ); /// 寻找轮廓 findContours( threshold_output, contours, hierarchy, CV_RETR_TREE, CV_CHAIN_APPROX_SIMPLE, Point(0, 0) ); /// 对每个轮廓计算其凸包 vector<vector<Point> >hull( contours.size() ); for( int i = 0; i < contours.size(); i++ ) {  
         convexHull( Mat(contours[i]), hull[i], false );//凸包计算
         }
    
       /// 绘出轮廓及其凸包
       Mat drawing = Mat::zeros( threshold_output.size(), CV_8UC3 );
       for( int i = 0; i< contours.size(); i++ )
          {
            Scalar color = Scalar( rng.uniform(0, 255), rng.uniform(0,255), rng.uniform(0,255) );
            drawContours( drawing, contours, i, color, 1, 8, vector<Vec4i>(), 0, Point() );//画轮廓
            drawContours( drawing, hull, i, color, 1, 8, vector<Vec4i>(), 0, Point() );//画凸包
          }
    
       /// 把结果显示在窗体
       namedWindow( "Hull demo", CV_WINDOW_AUTOSIZE );
       imshow( "Hull demo", drawing );
    凸包函数原型:void convexHull(InputArray points, OutputArray hull, bool clockwise=false, bool returnPoints=true);

    参数列表:2d点集合、输出的凸包、定向标识、操作标识;

    第一个参数:输入的2D点集合,存储在vector或者Mat中。

    第二个参数:输出的凸包,是一个包含索引的整数向量或者是点向量。在前者中,这个hull元素是在原始数组中的凸包点上基于0索引的;后者中,hull元素自身就是凸包点。(也就是说该参数要不是索引到原始图像找凸包点,或者是直接提取凸包点);

    第三个参数:定向标识,如果为true,那么输出的凸包就是顺时针方向的;不然就是逆时针方向的。假设的坐标系的x 轴指向右边,y 轴指向上面;

    第四个参数:操作标识,在矩阵Mat的情况中,当这个参数为true,该函数返回凸包点;否则返回指向凸包点的索引值;当在第二个参数是数组vector的情况中,该标识被忽略,输出是依赖vector的类型指定的,vector<int> 暗示returnPoints= true ,vector<Point>暗示 returnPoints= false。

        该函数使用Sklansky的算法来查找一个2d点集的凸包,在当前的执行情况下的复杂度是O(NlogN),可以参见convexhull.cpp中验证不同的函数变量的结果。

    3、给轮廓加上矩形或者圆形边界框


    这个还是挺有用的,有时候识别一个物体,想先用简单的数字图像处理方法得到更少的区域然后可以提取从而接着使用模式识别的方法进行训练和建立分类器。

    使用OpenCV函数 boundingRect 来计算包围轮廓的矩形框.

    使用OpenCV函数 minEnclosingCircle 来计算完全包围已有轮廓最小圆.

    int thresh = 100;
    int max_thresh = 255;
    RNG rng(12345);
    /// 转化成灰度图像并进行平滑 用来减少噪音点
      cvtColor( src, src_gray, CV_BGR2GRAY );
      blur( src_gray, src_gray, Size(3,3) );
    Mat threshold_output;
      vector<vector<Point> > contours;//存储轮廓点
      vector<Vec4i> hierarchy;//构建不同轮廓的层次结构
    
      /// 使用Threshold检测边缘
      threshold( src_gray, threshold_output, thresh, 255, THRESH_BINARY );
      /// 找到轮廓
      findContours( threshold_output, contours, hierarchy, CV_RETR_TREE, CV_CHAIN_APPROX_SIMPLE, Point(0, 0) );
    
      /// 多边形逼近轮廓 + 获取矩形和圆形边界框
      vector<vector<Point> > contours_poly( contours.size() );//创建contours.size()个空的多边形
      vector<Rect> boundRect( contours.size() );//创建contours.size()个矩形框
      vector<Point2f>center( contours.size() );//创建contours.size()个圆心
      vector<float>radius( contours.size() );//创建contours.size()个半径
    
      for( int i = 0; i < contours.size(); i++ )
         { 
           approxPolyDP( Mat(contours[i]), contours_poly[i], 3, true );//多边形逼近
           boundRect[i] = boundingRect( Mat(contours_poly[i]) );//获取某个轮廓的矩形框
           minEnclosingCircle( contours_poly[i], center[i], radius[i] );//生成某个轮廓的最小包含圆的圆心和半径
         }
    
    
      /// 画多边形轮廓 + 包围的矩形框 + 圆形框
      Mat drawing = Mat::zeros( threshold_output.size(), CV_8UC3 );
      for( int i = 0; i< contours.size(); i++ )
         {
           Scalar color = Scalar( rng.uniform(0, 255), rng.uniform(0,255), rng.uniform(0,255) );//随机颜色
           drawContours( drawing, contours_poly, i, color, 1, 8, vector<Vec4i>(), 0, Point() );//画轮廓-多边形
           rectangle( drawing, boundRect[i].tl(), boundRect[i].br(), color, 2, 8, 0 );//画矩形
           circle( drawing, center[i], (int)radius[i], color, 2, 8, 0 );//画圆形
         }
    
      /// 显示在一个窗口
      namedWindow( "Contours", CV_WINDOW_AUTOSIZE );
      imshow( "Contours", drawing );
    多边形逼近的函数原型:void approxPolyDP(InputArray curve, OutputArray approxCurve, double epsilon, bool closed);

    参数列表:输入数组、逼近的结果、逼近的精度、是否闭合的标识;

    第一个参数:一个2d点的输入向量可以存储在:

                              – std::vector or Mat    (C++ interface)
                             – Nx2 numpy array      (Python interface)
                            – CvSeq or ‘‘ CvMat      (C interface)

    第二个参数:逼近的结果,该参数的类型应该匹配输入曲线的类型。(在c接口中,逼近的曲线存储在内存存储区中,返回的是指向该内存的指针,我列出来的都是cpp接口,所以该句可忽略);

    第三个参数:指定逼近精度的参数。这是介于原始曲线与它的逼近之间的最大的距离;

    第四个参数:如果为真,逼近的曲线是闭合的(即,将首尾连起来);否则就不是闭合的。

    该函数使用一个更少顶点的曲线/多边形来逼近另一个曲线或多边形,使得介于他们之间的距离小于或者等于指定的精度。该函数使用的是Douglas-Peucker算法。http://en.wikipedia.org/wiki/Ramer-Douglas-Peucker_algorithm


    计算一个点集合的up-right边界矩形的函数原型:Rect boundingRect(InputArray points);

    参数列表:输入点;

    第一个参数:输入的2D点集合,存储在vector或者Mat中。

    该函数对指定的点集合计算并返回最小的up-right矩形。


    计算一个2d点集合的最小区域圆的函数原型:void minEnclosingCircle(InputArray points, Point2f& center, float& radius)

    参数列表:输入的点集合、输出的圆心坐标、输出的半径;

    第一个参数:2D点的输入向量,存储在:

                            – std::vector<> or Mat (C++ interface)
                             – CvSeq* or CvMat* (C interface)
                            – Nx2 numpy array (Python interface)

    第二个参数:圆的输出的圆心;

    第三个参数:圆的输出的半径

    该函数在一个2D点集合中使用一个迭代算法来查找最接近的圆,见minarea.cpp中有例子。

    4、给轮廓加上倾斜的边界框和椭圆

    该部分与上面部分很相似,不过采用的函数是不同的。


    int thresh = 100;
    int max_thresh = 255;
    RNG rng(12345);
     /// 转为灰度图并模糊化 来减小噪音点
      cvtColor( src, src_gray, CV_BGR2GRAY );
      blur( src_gray, src_gray, Size(3,3) );
    
     Mat threshold_output;
      vector<vector<Point> > contours;//存储轮廓的向量
      vector<Vec4i> hierarchy;//轮廓的层次结构
    
      /// 阈值化检测边界
      threshold( src_gray, threshold_output, thresh, 255, THRESH_BINARY );
      /// 寻找轮廓
      findContours( threshold_output, contours, hierarchy, CV_RETR_TREE, CV_CHAIN_APPROX_SIMPLE, Point(0, 0) );
    
      /// 对每个找到的轮廓创建可倾斜的边界框和椭圆
      vector<RotatedRect> minRect( contours.size() );//存储contours.size()个旋转的边界框***************此处是重点,RotateRect
      vector<RotatedRect> minEllipse( contours.size() );////存储contours.size()个旋转的椭圆*************
    
      for( int i = 0; i < contours.size(); i++ )
         { 
            minRect[i] = minAreaRect( Mat(contours[i]) );//获取包含最小区域的矩形
           if( contours[i].size() > 5 )
             { 
             minEllipse[i] = fitEllipse( Mat(contours[i]) );//获取椭圆所需的信息
             }
         }
    
      /// 绘出轮廓及其可倾斜的边界框和边界椭圆
      Mat drawing = Mat::zeros( threshold_output.size(), CV_8UC3 );
      for( int i = 0; i< contours.size(); i++ )
         {
           Scalar color = Scalar( rng.uniform(0, 255), rng.uniform(0,255), rng.uniform(0,255) );
           // contour
           drawContours( drawing, contours, i, color, 1, 8, vector<Vec4i>(), 0, Point() );//画对象的轮廓
           // ellipse
           ellipse( drawing, minEllipse[i], color, 2, 8 );//画椭圆
           // rotated rectangle
           Point2f rect_points[4]; 
           minRect[i].points( rect_points );//将旋转矩形的顶点赋值给实参;
           for( int j = 0; j < 4; j++ )
              line( drawing, rect_points[j], rect_points[(j+1)%4], color, 1, 8 );//画出倾斜矩形,采用的是画四条线的形式实现的
         }
    
      /// 结果在窗体中显示
      namedWindow( "Contours", CV_WINDOW_AUTOSIZE );
      imshow( "Contours", drawing );
    在输入的2D点集合中找到一个最小区域的旋转后的矩形的函数原型:RotatedRect minAreaRect(InputArray points)
    参数列表:输入的2D点集合
    第一个参数:输入的2D点的向量,存储在:
                                   – std::vector<> or Mat (C++ interface)

                                  – CvSeq* or CvMat* (C interface)
                                  – Nx2 numpy array (Python interface)

    该函数计算并返回一个指定的点集合的最小区域边界的矩形,在minarea.cpp中有例子,开发者应该注意当数据是接近包含的矩阵元素边界的时候 返回的旋转矩形可以包含负指数。


    沿着一个2D点集合进行拟合一个椭圆的函数原型:RotatedRect fitEllipse(InputArray points)

    第一个参数:输入的2D点集合,可以存储在:

                          – std::vector<> or Mat (C++ interface)
                         – CvSeq* or CvMat* (C interface)
                         – Nx2 numpy array (Python interface)

    该函数计算在对一个2D点集合拟合的最好的椭圆(最小二乘作为损失函数来判别)。它返回旋转的矩形中内含的椭圆,使用的算法是[FItzgibbon95].开发者应该注意当数据点很靠近包含的矩阵元素的边界的时候,返回的椭圆/旋转矩形数据包含负指数例子代码在fitellipse.cpp中。

    5、计算轮廓的矩

    使用OpenCV函数 moments 计算图像所有的矩(最高到3阶)

    使用OpenCV函数 contourArea 来计算轮廓面积

    使用OpenCV函数 arcLength 来计算轮廓或曲线长度

    int thresh = 100;
    int max_thresh = 255;
    RNG rng(12345);
     /// 把原图像转化成灰度图像并进行平滑
      cvtColor( src, src_gray, CV_BGR2GRAY );
      blur( src_gray, src_gray, Size(3,3) );
    Mat canny_output;
      vector<vector<Point> > contours;//存储轮廓
      vector<Vec4i> hierarchy;//存储轮廓的层次结构
    
      /// 使用Canndy检测边缘
      Canny( src_gray, canny_output, thresh, thresh*2, 3 );
      /// 找到轮廓
      findContours( canny_output, contours, hierarchy, CV_RETR_TREE, CV_CHAIN_APPROX_SIMPLE, Point(0, 0) );
    
      /// 计算矩
      vector<Moments> mu(contours.size() );
      for( int i = 0; i < contours.size(); i++ )
         { 
             mu[i] = moments( contours[i], false );//计算轮廓的矩
          }
    
      ///  计算中心矩:
      vector<Point2f> mc( contours.size() );
      for( int i = 0; i < contours.size(); i++ )
         { 
              mc[i] = Point2f( mu[i].m10/mu[i].m00 , mu[i].m01/mu[i].m00 ); 
          }
    
      /// 绘制轮廓
      Mat drawing = Mat::zeros( canny_output.size(), CV_8UC3 );
      for( int i = 0; i< contours.size(); i++ )
         {
           Scalar color = Scalar( rng.uniform(0, 255), rng.uniform(0,255), rng.uniform(0,255) );
           drawContours( drawing, contours, i, color, 2, 8, hierarchy, 0, Point() );
           circle( drawing, mc[i], 4, color, -1, 8, 0 );
         }
    
      /// 显示到窗口中
      namedWindow( "Contours", CV_WINDOW_AUTOSIZE );
      imshow( "Contours", drawing );
    
      /// 通过m00计算轮廓面积并且和OpenCV函数比较
      printf("\t Info: Area and Contour Length \n");
      for( int i = 0; i< contours.size(); i++ )
         {
           printf(" * Contour[%d] - Area (M_00) = %.2f - Area OpenCV: %.2f - Length: %.2f \n", 
                  i, mu[i].m00, contourArea(contours[i]), arcLength( contours[i], true ) //第几个轮廓,轮廓的空间矩,轮廓的面积,轮廓的周长
                   );//计算轮廓面积和 轮廓或曲线的长度
           Scalar color = Scalar( rng.uniform(0, 255), rng.uniform(0,255), rng.uniform(0,255) );
           drawContours( drawing, contours, i, color, 2, 8, hierarchy, 0, Point() );
           circle( drawing, mc[i], 4, color, -1, 8, 0 );
     

    计算一个多边形或者光栅形状所有的矩最高到3阶的函数原型:Moments moments(InputArray array, bool binaryImage=false )

    参数列表:输入数组、标识

    第一个参数:光栅图像(单通道,8-bit或者浮点2D数组);或者一个2D点(Point 或者 Point2f)的数组(1xN或者Nx1);

    第二个参数:如果为真,所有的非0图像像素被视为 1 ;该参数只用在图像上。

    该函数计算一个向量形状或者光栅形状的矩到3阶。结果返回在结构体Moments。

    class Moments
    {
    public:
    <span style="white-space:pre">		</span>Moments();
    <span style="white-space:pre">		</span>Moments(double m00, double m10, double m01, double m20, double m11,
    <span style="white-space:pre">		</span>double m02, double m30, double m21, double m12, double m03 );
    <span style="white-space:pre">		</span>Moments( const CvMoments& moments );
    <span style="white-space:pre">		</span>operator CvMoments() const;
    // spatial moments
    <span style="white-space:pre">		</span>double m00, m10, m01, m20, m11, m02, m30, m21, m12, m03;
    // central moments
    <span style="white-space:pre">		</span>double mu20, mu11, mu02, mu30, mu21, mu12, mu03;
    // central normalized moments
    <span style="white-space:pre">		</span>double nu20, nu11, nu02, nu30, nu21, nu12, nu03;
    }
    在光栅图像的时候,空间矩以下面的方式计算:

     

    中心矩以下面的方式计算:

    这里是团块中心:

    归一化的中心矩按照下面的方式计算:

    notes:
    轮廓的矩以同样的方式定义,不过是使用Green(http://en.wikipedia.org/wiki/Green_theorem).)的公式计算的。所以因为一个受限的光栅分辨率,轮廓的矩的计算轻微的不同于同样的光栅轮廓的矩计算。

    因为轮廓矩是使用Green公式计算的,所以可能看到奇怪的与自相交的轮廓结果,例如:一个蝴蝶形状轮廓的0区域(m00)。


    计算一个轮廓的面积的函数原型:double contourArea(InputArray contour, bool oriented=false )

    第一个参数:输入的2D点的向量(轮廓索引),存储在vector或者Mat中;

    第二个参数:定向区域标识。如果为true,该函数返回一个有符号的面积值,依赖于轮廓的方向(顺时针或者逆时针)。使用这个特征,你可以通过使用一个区域的符号决定一个轮廓的方向。默认情况下,该参数是false的,这意味着会返回绝对值。

    该函数计算一个轮廓的面积,相似于moment()。该面积是使用Green公式计算的。所以,如果你通过使用drawContours()或者fillPoly()画这个轮廓,返回的面积和非0像素的个数,会不一样。同样的该函数差不多会在一个有着自相交的轮廓上给出一个错误的结果。(也就是这时候该函数不靠谱)。

    例子代码:

       vector<Point> contour;
    contour.push_back(Point2f(0, 0));
    contour.push_back(Point2f(10, 0));
    contour.push_back(Point2f(10, 10));
    contour.push_back(Point2f(5, 4));
       double area0 = contourArea(contour);
    vector<Point> approx;
       approxPolyDP(contour, approx, 5, true);
       double area1 = contourArea(approx);
    cout << "area0 =" << area0 << endl <<
    <span style="white-space:pre">	</span>"area1 =" << area1 << endl <<
    <span style="white-space:pre">	</span>"approx poly vertices" << approx.size() << endl;

    计算一个轮廓的周边或者一个曲线长度的函数原型:double arcLength(InputArray curve, bool closed)

    参数列表:

    第一个参数:输入的2D点向量,存储在vector或者Mat中;

    第二个参数:表示该曲线是否闭合。

    该函数计算一个曲线的长度或者一个闭合轮廓的周长。

    6、多边形测试


    是执行一个轮廓中点的测试。

     /// 创建一个图形     const int r = 100;
      Mat src = Mat::zeros( Size( 4*r, 4*r ), CV_8UC1 );//创建画布
    
      /// 绘制一系列点创建一个轮廓:
      vector<Point2f> vert(6);//一个有6个点的轮廓
    
      vert[0] = Point( 1.5*r, 1.34*r );
      vert[1] = Point( 1*r, 2*r );
      vert[2] = Point( 1.5*r, 2.866*r );
      vert[3] = Point( 2.5*r, 2.866*r );
      vert[4] = Point( 3*r, 2*r );
      vert[5] = Point( 2.5*r, 1.34*r );//自己定轮廓的6各点
    
      /// 在src内部绘制
      for( int j = 0; j < 6; j++ ){ 
          line( src, vert[j],  vert[(j+1)%6], Scalar( 255 ), 3, 8 ); //将该轮廓画出来
          }
    
      /// 得到轮廓
      vector<vector<Point> > contours; vector<Vec4i> hierarchy;
      Mat src_copy = src.clone();
    
      findContours( src_copy, contours, hierarchy, RETR_TREE, CHAIN_APPROX_SIMPLE);//查找轮廓
    
      /// 计算到轮廓的距离
      Mat raw_dist( src.size(), CV_32FC1 );
    
      for( int j = 0; j < src.rows; j++ ){ 
        for( int i = 0; i < src.cols; i++ ){ 
               raw_dist.at<float>(j,i) = pointPolygonTest( contours[0], Point2f(i,j), true ); 整个画布,逐点测试到第0个轮廓的距离
              }
         }
    
      double minVal; double maxVal;
      minMaxLoc( raw_dist, &minVal, &maxVal, 0, 0, Mat() );//查找全局最大最小
      minVal = abs(minVal); maxVal = abs(maxVal);
    
      /// 图形化的显示距离
      Mat drawing = Mat::zeros( src.size(), CV_8UC3 );
    
      for( int j = 0; j < src.rows; j++ ){ 
          for( int i = 0; i < src.cols; i++ ) {
                if( raw_dist.at<float>(j,i) < 0 ){ 
                   drawing.at<Vec3b>(j,i)[0] = 255 - (int) abs(raw_dist.at<float>(j,i))*255/minVal;
                 }
                else if( raw_dist.at<float>(j,i) > 0 ){
                   drawing.at<Vec3b>(j,i)[2] = 255 - (int) raw_dist.at<float>(j,i)*255/maxVal;
                  }
                else{ 
                  drawing.at<Vec3b>(j,i)[0] = 255; drawing.at<Vec3b>(j,i)[1] = 255; drawing.at<Vec3b>(j,i)[2] = 255;
                  }
              }
         }
    
      /// 创建窗口显示结果
      char* source_window = "Source";
      namedWindow( source_window, CV_WINDOW_AUTOSIZE );
      imshow( source_window, src );
      namedWindow( "Distance", CV_WINDOW_AUTOSIZE );
      imshow( "Distance", drawing );
    
    执行一个轮廓中点的测试的函数原型:double pointPolygonTest(InputArray contour, Point2f pt, bool measureDist)
    参数列表:输入的轮廓、测试的点、标识;

    第一个参数:输入的轮廓;

    第二个参数:针对某个的测试的点;

    第三个参数:如果为真,该函数将评估一个点到最近的轮廓的边界的有符号距离;否则只检查该点是否在轮廓内部。

    该函数决策着该点在一个轮廓的内部,外部,或者刚好在轮廓的边界上(或者刚好在轮廓的顶点上)。该函数返回为正(内部),负(外部),0(在边界上)。当measureDist = false,该函数返回+1,-1,0;不然返回的值是一个点到最近轮廓边的有符号距离值。


    图像处理部分的例子算是都过了一遍了,接下来是ml部分

    展开全文
  • 在上一篇博客(OCR识别字符排列圆形或字体倾斜处理办法)中我们分析了如何矫正倾斜的字符,这里直接上代码 * 加载图片,注意更改文件名 read_image (Image, 'ImageName') text_line_slant (Image, Image, 60, ...
  • 分类: 图像处理算法2009-06-20 20:5454441人阅读评论(10)收藏举报 算法blog活动 识别算法概述: SIFT/SURF基于灰度图, 一、首先建立图像金字塔,形成三维的图像空间,通过Hessian矩阵获取每一层的局部极大值,然后...
  • 文本图像畸变矫正

    2020-06-01 14:47:01
    二十世纪六十年代兴起的OCR(Optical Character Recognition,光学字符识别)技术,使得文档能以图像的形式被分析与识别,一定程度上实现了文本识别的自动化。然而,文档图像识别效果的优劣与其质量有着密切的联系,...
  • 文章目录一、文本检测概述二、Pixel-Anchor 网络详解2.1、Pixel-Anchor网络结构2.2、像素级别语义分割模块(Pixel based Module)2.3、锚检测回归模块(Anchor based Module)2.4、后处理2.5、Pixel-Anchor检测效果...
  • OCR技术概览

    2019-09-11 15:59:57
    , 印刷体识别比手写体识别要简单, 因为印刷体更规范, 字体来自于计算机字库, 尽管印刷过程中可能会发生不清晰\粘连, 这些都可以通过一些"腐蚀"/"膨胀"图像处理技术还原, 但是手写体由于个体差异的存在, 还是非常难的...
  • 一个轮廓一般对应一系列的点,也就是图像中的一条曲线.表示的方法可能根据不同情况而有所不同.有多重方法可以表示曲线.在OpenCV中一般用序列来存储轮廓信息.序列中的每一个元素是曲线中一个点的位置.关于序列表示的...
  • pinhole camera 1st person camera world 定义了第1人称相机世界的坐标系 定义第1人称相机世界 我称它第1人称相机世界是...然后想象y轴是竖直向下的 如果使用右手法则 将x叉乘y 它是指向下的 你的拇指将指向z方向...
  • 【珍藏孤版】摄影技能技巧大全   (2013-07-13 12:19:29) 转载▼ 标签:  转载   学而时习之,不亦说乎. 原文地址:【珍藏孤版】摄影技能技巧...光圈配快门,曝光要先练,
  • (1)首先我们来看车牌部分与其他部分的区别,及车牌的特征,整体的车牌定位及提取方案即时基于此(基于灰度图像中考虑): a.车牌部分是矩形 b.车牌具有特定的长宽比; c.车牌的面积一定 (2)接着我们对整幅...
  • 人类是如何识别一个物体的呢,当然要对面前的这个物体为何物要有一个概念,人类一生下来就开始通过视觉获取世间万物的信息,包括一种物体形状、颜色、成分等,以及通过学习认识到这种物体的其他信息比如物理的、化学...
  • 有兴趣的读者可以从本书开始,通过图像分类、识别、检测和分割的案例,逐步深入卷积神经网络的核心,掌握深度学习的方法和精髓,领会 AlphaGo 战胜人类世界冠军的奥秘。 作者简介 李玉鑑(鉴) 北京工业大学教授,...
  • 现在只是得到MPU6050的一些原始数据,还未做滤波处理。 接下来先讲,加速度计和陀螺仪的计算公式,然后进一步延伸出姿态滤波。 一、加速度计 (1)计算公式 参看:Arduino教程:MPU6050的数据获取、分...
1 2 3 4 5 ... 10
收藏数 187
精华内容 74