2013全国数学建模竞赛 ——碎纸片拼接复原.doc

2013高教社杯全国大学生数学建模竞赛重庆工商大学 姜木北小组作品编 号 专 用 页赛区评阅编号（由赛区组委会评阅前进行编号）：赛区评阅记录（可供赛区评阅时使用）：评阅人评分备注全国统一编号（由赛区组委会送交全国前编号）：全国评阅编号（由全国组委会评阅前进行编号）碎纸片的拼接复原摘 要目前，“碎片拼接复原”技术在司法物证复原、历史文物修复及社会生活各项领域扮演着重要角色，对于碎片数量特别巨大而人工又难以在短时间内完成碎片拼接时，要找到一种高效快捷的自动拼接方法已变得尤为重要。本文针对只有中英文的碎片拼接问题，综合分析了从单一的纵切到纵横切以及纵横切双面碎片这三个不同的情况，提出了碎片拼接复原的解决方案.在问题一中，对于仅有“纵切”且数量相对较少的碎纸片，我们基于边缘去噪和采用构建碎纸图片的左右边缘二值矩阵提取相似度分析的方法，再通过两张图片左右相似度匹配排序，得到附件1和附件2中的碎纸排序（见表2和表3），并运用Matlab的图像处理工具箱，按排列顺序导入碎纸片得到相应拼接结果（见附录附件一）.在问题二中，由于碎纸片数量相对较多，同时存在横切和纵切的情况，在问题一的基础上增加了碎纸片的上下边缘相似度匹配。在进行人工干预，找到第一张起始碎纸片作为匹配起点后，我们基于索贝尔算子的原理，对碎纸片灰度值进行边缘相似度的旋转检测和比较匹配，最后进行二叉树搜索排序（见表4和表5）。对附件3和4的碎纸图片拼接出的结果详见附录中的附件二.在问题三中，由于碎纸片是两面的并且碎纸片数量更多，若采用第二问的求解方案则加大了求解难度同时也存在较大误差。因此，我们基于蚁群算法（ACA）的SIFT特征点匹配原理来求解。先提取碎纸图片特征点，然后基于蚁群算法的最优化快速比对匹配，最后基于ACA的搜索排序对碎纸片拼接。Matlab编程所求得的排序结果详见表6和表7，附件5中的碎纸片拼接复原结果见附录中的附件三.在问题的解决中，我们得出结论:碎纸图片导入量越小，图片匹配出的效果越佳，在相似度的匹配上，出现的误差减小,最后拼出的图像效果好，人工干预量也相对小。本建模考虑到了图片噪声对图片拼接时的影响，选择了去噪效果较好的邻域平均法对图片进行处理。但是，为了解题方便，我们忽略了碎纸机切纸时碎纸片可能产生的边界遗失破损。在这里我们的改进是，对于边界遗失图像碎片的修复，根据复原后的整体形状，可以根据线连续性来拟合此类线段,从而得到较为完整的图像.关键词：碎纸拼接、腐蚀算法、蚁群算法、图像特征匹配、邻域平均法一、问题的重述破碎文件的拼接在司法物证复原、历史文献修复以及军事情报获取等领域都有着重要的应用。传统上，拼接复原工作需由人工完成，准确率较高，但效率很低。特别是当碎片数量巨大，人工拼接很难在短时间内完成任务。随着计算机技术的发展，人们试图开发碎纸片的自动拼接技术，以提高拼接复原效率。现解决如下问题：1. 对于给定的来自同一页印刷文字文件的碎纸机破碎纸片（仅纵切），建立碎纸片拼接复原模型和算法，并针对附件1、附件2给出的中、英文各一页文件的碎片数据进行拼接复原。如果复原过程需要人工干预，请写出干预方式及干预的时间节点。复原结果以图片形式及表格形式表达.2. 对于碎纸机既纵切又横切的情形，请设计碎纸片拼接复原模型和算法，并针对附件3、附件4给出的中、英文各一页文件的碎片数据进行拼接复原。如果复原过程需要人工干预，请写出干预方式及干预的时间节点。复原结果表达要求同上.3. 上述所给碎片数据均为单面打印文件，从现实情形出发，还可能有双面打印文件的碎纸片拼接复原问题需要解决。附件5给出的是一页英文印刷文字双面打印文件的碎片数据。请尝试设计相应的碎纸片拼接复原模型与算法，并就附件5的碎片数据给出拼接复原结果，结果表达要求同上.二、问题的分析本题是一个关于碎纸片拼接复原的问题，针对题目的要求，我们对于“碎纸片拼接复原”等一系列问题的处理特点和处理方案，就每个问题做除了以下的具体分析： 21针对问题1的分析结合给出的碎纸机碎片情况，分别是一页中文和英文文件的碎片数据，也就是说，在每张碎片中都会出现相当多的完整或不完整的汉字（英文）或者结构，要解决把纵切的碎纸片拼接复原的问题：首先，利用MATLAB编程建立图片的二值矩阵。鉴于纸片的大小相同，并且只有黑白两色，通过观察附件一的碎片，每张碎片最左端或者最右端不是全空，这些碎片就可以全部转化为二值图，0代表白色，1代表黑色。其次，采用邻域平均法对图像进行边缘去噪，进一步利用分割最佳阈值的迭带算法来达到消除或减少噪声影响的效果。最后，比较图像左右边缘的相似情况,从而得出两个边缘图像的相似度,以此来度量两个图像能否匹配。复原过程中需要人工干预，将依顺序显示的未拼接碎片，采用半自动拼接，即人工选择碎片并拼接到计算机屏幕上。将序号排列好后；再通过基于MATLAB编程将19张图片分别排好序通过调用imshow()函数合成图片，这样就能得到最终求解结果.22针对问题2的分析求解要求中我们可以知道该问题在第一问的基础上增加了“横切”碎片，出现的纸张既有横切又有纵切。解决方案与第一问有所不同；增加了“横切”过后，在第一问的左右相似度比较之上增加了上下边缘的相似度比较和匹配。但考虑到求解精度我们采用灰度值（0255）求解的图片处理方法；同样基于MATLAB编程对附件3和附件4中的碎片分别进行灰度处理得到图片矩阵的灰度值；并采用基于Sobel 的算法对个碎片进行四个边的边缘化处理，提取边缘特征，并对其进行边缘相似度的计算。同时对边缘化处理后的灰度值进行聚类分析确定一个阀值，也就是确定一个误差值。只考虑切缝处并对其相似度一致的缝合。又由于该问中需要拼接的图片数量众多，我们采用基于二叉树的搜索排序算法对相似度矩阵排序并得到该相似度矩阵所对应的图片序列号数组，该排列的序列号数组即为碎片的组合方式；最后采用MATLAB循环调用序列号数组对应的碎片名拼接出所求结果达到解决该问题的目的.2.3针对问题3的分析对于问题3其实是在问题二的基础上在增加了正反面的拼接。另外通过分析题目中告诉了的类容可以发现附件中同一个序号表示是一张碎纸图片，a和b表示同一张碎片的了两面，但并不确定a和b的正反。这样同第二问中的解决思路相似，但由于问题三中对应的附件5的碎图片数据量过大如若采用第二问的求解方式也并不理想所以我们采用基于蚁群算法（ACA）的SIFT特征点图像碎片拼接方法。针对附件五中的图片基于MATLAB采用所设计的算法提取其特征点并建立特征点匹配矩阵。再利用这些特征点匹配矩阵根据理论建模设计基于MATLAB的算法求解出顺序排列的部分序列号和误差节点。将其中求得的序列号对应的碎纸图片导入MATLAB中求解出拼接图像，最后通过人工干预的方式补全各个残缺节点的碎片得到完整的纸片图像已达到拼接出碎纸片的正反两面的目的.三、模型假设1、 假设未碎纸张的文字行方向沿水平方向，字与字之间有间隔且字宽度与高度比值1/3；2、 假设碎片模型为理想模型，碎片表面光滑平整无磨损且厚度为零；3、 假设在切割过程中除边界外，其他切割线都切割文字；4、 假设每张碎纸片的大小一致；5、 假设碎纸机切的纸片无损坏；6、 假设碎图片中的文字符字体格式相同，英文字符字体格式相同. 四、符号的定义图像平面上的一系列平均噪声点，第个二值矩阵图像的分割阀值，即灰度值阀值，灰度均值，待测边缘图像边缘点经横（纵向）边缘检测的图像梯度方向第i像素和它的K连通像素边缘相似度不同中心像素点的协议变量，为启发式因子的相对重要程度为信息素的相对重要程度为信息素蒸发系数表示信息素的持久性系数为窗口信息素含量五、模型的建立及求解5.1 问题一模型的建立及求解算法的流程图：载入图像建立图片的二值矩阵采用邻域平均法进行边缘去噪比较图像左右边缘的相似情况运用腐蚀算法轮廓提取方法来消除图像的左右边界的作用图15.1.1、把图像转化为二值矩阵二值图像是指每个像素不是黑就是白，其灰度值没有中间过渡的图像。二值图像一般用来描述文字或者图形，其优点是占用空间少，人们经常用黑白、B&W、单色图像表示二值图像，但是也可以用来表示每个像素只有一个采样值的任何图像，例如灰度图像等.因此针对问题一，我们先将对应的碎纸图片用Matlab调用im2bw函数对导入的碎纸图片求解其二值矩阵.5.1.2、图像边缘去噪邻域平均法邻域平均法是一种利用Box模板对图像进行模板操作（卷积操作）的图像平滑方法最简单的平滑滤波是将原图中一个像素的灰度值和它周围邻近像素的灰度值相加，然后将求得的平均值作为新图中该像素的灰度值。它采用模板计算的思想，模板操作实现了一种邻域运算，即某个像素点的结果不仅与本像素灰度有关，而且与其邻域点的像素值有关.将附件1或2中已经进行二值化的碎纸片矩阵惊醒边缘去噪处理，消弱噪声点起到平滑作用，减小相似度误差，增大匹配准确度.设为给定的含有噪声的图像，经过邻域平均处理后的图像为，则邻域平均法也可以用数学公式表达： (1)式中：，；是以为中心的邻域的集合，M是S内的点数.邻域平均法的思想是通过一点和邻域内的像素点求平均来去除突变的像素点，从而滤掉一定的噪声，其主要在实际应用中，也可以根据不同的需要选择使用不同的模板尺寸，如33、55、77、99等。邻域平均处理方法是以图像模糊为代价来减小噪声的，且模板尺寸越大，噪声减小的效果越显著。如果是噪声点，其邻近像素灰度与之相差很大，采用邻域平均法就是用邻近像素的平均值来代替它，这样能明显消弱噪声点，使邻域中灰度接近均匀，起到平滑灰度的作用.从理论上看: a、要改变均值滤波的滤波能力，就要改变模板的大小，模板越大去噪能力越强，反之越弱； b、由于是在整幅图上以取局部灰度平均的方法来淡化噪声，所以不可避免的会使原有图像的清晰度、对比度下降，会损失原有图像的细节，模板越大损失的越多； c、对整幅图像用同一个模板，不能适应局部具体情况，会出现有滤波不均匀的情况.5.1.3、分割阀值的迭代分析法 阀值分割方法是把图像的灰度分成不同的等级，然后用设置灰度门限的方法确定欲分割物体的边界。当用阀值来分割口标与背景时，如果某一灰度值g是某图像的分割阀值，即小于g的灰度点将构成口标(不妨如此假设)，而大于g的灰度点就构成背景一般而言，如果灰度值g可以作为图像的一个阀值，那么它应该使按这个阀值划分口标和背景的错误分割的图像像素点数为最小.如果前景物体的内部具有均匀一致的灰度值，并分布在另一个灰度值的均匀背景上，那么图像的灰度直方图应具有明显的双峰。可是在许多情况下，噪声的干扰使峰谷的位置难以判定或者结果不稳定本文采用迭代算法，有效地消除或减少噪声对灰度门限值g的影响.设有一幅混入噪声的图像g (x),力是由原始图像f (x,y)和e(x,y)叠加而成的。即: （2）这里假设各点的噪声是互不相关的，且具有零均值，标准差为。通过阀值分割将图像分割为两部分，由于噪声是随机作用于图像的像素点上，则可以认为在分割出口标和背景图像上噪声干扰仍为，即: （3） （4）在迭代算法中，需要对分割出的图像分别求其灰度均值，则: （5） （6）上式说明，随着迭代次数的增加，平均灰度值将趋向于真值。因此，用迭代算法求得的最佳阀值不受噪声干扰的影响根据上述分析，针对该问题的迭代算法描述如下.（1）首先选择一个近似阀值作为估训一值的初始值，然后进行分割，产生子图像，并根据子图像的特性来选择新的阀值，再用新的阀值分割图像，经过几次循环，使错误分割的图像像素点降到最少。这样做的效果好于用初始阀值直接分割图像的效果，阀值的改进策略是迭代算法的关键。算法步骤如下: 选择一个初始阀值的估算值 （7）式中，分别表示图像中的最小和最大灰度值.（2）利用阀值把图像分割成两组，和，其中： （8） （9）(3) 计算区域和的灰度值均值，其中： （10） （11）式中，f (x,y)是图像上(i，j)点的灰度值，N(i,j)是(i,j)点的权重系数，一般N(i,j)=1.0。（4）选择新的阈值 （12）（5）如果，则结束,否则,转步骤(2)。 经过图像分割处理后,有效减少了噪声对图像的干扰范围,并使图像的边缘邻域像素点离噪声干扰的敏感区,从而提高了图像边缘检测抗干扰能力.5.1.4、腐蚀算法轮廓提取方法经过上述三求解阶段过后，图像轮廓提取算法就变得非常简单。基于二值图像轮廓提取的算法就是掏空内部点,如果原图中有一点为黑,且它的8个相领点都是黑色时,判定该点是原图像的内部点,则将该点删除。最后,经过这样的算法处缘检测4, 5。在数学形态学运算中,腐蚀具有消除物体边界点的作用。结构元素取33的黑点块,腐蚀将使物体的边界沿周边减少一个像素。那么边缘检测实际上相当于用33块的9个点结构元素对原图进行腐蚀,再用原图像减去腐蚀的图像。令X为图像,B为结构元素,BZ表示结构元素B平移Z后的结果, 代表结构元素关于原点的对称集合。其数学表示如下: （13）则腐蚀的运算定义: (14)式中, 的结构元素; ED(X):代表图像X的边界.该方法检测到的物体边缘宽度仅为一个像素,因此具有较高的定位精度。同时,二值化处理后的图像具有完整的轮廓,所以轮廓提取检测到的边缘具有连续性.5.1.5、边缘相似度 边缘相似度的基本思想是用标准的无噪图像产生出标准的边缘图像，设为，将含噪图像用噪声抑制算法滤波处理后产生待测边缘图像，设为,比较两个边缘图像上各个边缘点的相似情况，从而得出两个边缘图像的相似度，以此来度量滤波算法的边缘保持能力. 对于上的一个边缘点，其边缘相似度定义为： （15）由上式可见，单点边缘相似度反映的是检测边缘点与对应的原标准边缘点之间的位置偏差:如果完全重合，则为最大相似度1;如果相差k个像素点，则其下降为，如果是一个误检点包括错检点和漏检点，即原来不是边缘点，但检测成了一个边缘点；或者原来有边缘，但是没检测出来，则为-1.整幅图像的边缘相似度则定义为所有边缘点相似度的平均公式： （16）式中M为中的边缘点数.5.1.6、基于计算机MATLAB的模型求解（1）计算机编程求解得到图片的二值矩阵该问题所见模型基于计算机仿真的MATLAB求解碎纸片的二值矩阵可以得到一个的0、1矩阵列.以对附件一的纵切碎纸片拼接为例,罗列出前6张碎纸图片求得的二值矩阵如下：、表1: 二值矩阵000.bmp 001.bmp 002.bmp000011000110011100111101100101010000001011010100101000101000100010110111003.bmp 004.bmp 005.bmp010010110100111101011001000100101011010100101000101000001010011111010101（2）相似度匹配结果通过相似度匹配分析可以得到每张碎纸图片右边缘与另一张碎纸图片左边缘的相似度匹配分析。若高于90%则匹配吻合。以附件一为例：先人工干预确定第一张图片“008.bmp”为起始点，它的二值矩阵右边与“014.bmp”相似度匹配结果最高.（3）求得的排序结果表2：对于附件一基于相似度分析匹配求得的碎纸片自动排序8141215310216145913181171706表3：对于附件二基于相似度分析匹配求得的碎纸片自动排序4 7381619121621014119131518175（说明：“8”表示第九张图片“008.bmp”排在第一列；附件1和2碎纸图片的拼接结果详见附件）5.2 问题二模型建立及求解算法流程图：Sobel 的算法对个碎片进行四个边的边缘化处理对边缘化处理后的图像进行聚类分析确定阀值二叉树的搜索排序算法对图片排序二元相似性度量算法对图像进行分析导入图像图2在实际图像边缘检测问题中，图像的边缘作为图像的一种基本特征，经常被应用到较高层次的图像应用中去。图像边缘是图像最基本的特征之一，往往携带着一幅图像的大部分信息。而边缘存在于图像的不规则结构和不平稳现象中，也即存在于信号的突变点处，这些点给出了图像轮廓的位置，这些轮廓常常是我们在图像边缘检测时所需要的非常重要的一些特征条件，这就需要我们对一幅图像检测并提取出它的边缘.5.2.1、Sobel operator索贝尔算子（Sobel operator）是图像处理中的算子之一，主要用作边缘检测。在技术上，它是一离散性差分算子，用来运算图像亮度函数的梯度之近似值。在图像的任何一点使用此算子，将会产生对应的梯度矢量或是其法矢量.该算子包含两组3x3的矩阵，分别为横向及纵向，将之与图像作平面卷积，即可分别得出横向及纵向的亮度差分近似值。如果以A代表原始图像，Gx及Gy分别代表经横向及纵向边缘检测的图像，其公式如下: 图像的每一个像素的横向及纵向梯度近似值可用以下的公式结合，来计算梯度的大小。 (17)然后可用以下公式计算梯度方向 (18)在以上例子中，如果以上的角度等于零，即代表图像该处拥有纵向边缘，左方较右方暗。在边沿检测中，常用的一种模板是Sobel 算子。Sobel 算子有两个，一个是检测水平边沿的；另一个是检测垂直平边沿的。Sobel算子对于象素的位置的影响做了加权，因此效果更好。Sobel算子另一种形式是各向同性Sobel(Isotropic Sobel)算子，也有两个，一个是检测水平边沿的，另一个是检测垂直平边沿的。各向同性Sobel算子和普通Sobel算子相比，它的位置加权系数更为准确，在检测不同方向的边沿时梯度的幅度一致.5.2.2、调整簇阀值的加速聚类方法调整簇阀值的加速聚类方法的基本思想是:设置簇调整阀值,在每次迭代完成后,只有簇调整量大于阀值的簇才参与下一次迭代调整。簇调整阀值可以设置为簇中心迁移的调整阀值和簇中样本数目调整阀值二种情况:(1)记某一簇在第1次迭代后聚类中心与初始的簇中心之间的距离为,在第次迭代后新的簇中心与原簇中心之间的距离,如果,则以为代表的簇不再参与下轮迭代调整,否则,继续参与迭代调整.(2)对于任一个簇,假设样本数目,在每次迭代完成后,该簇中样本调整个数为(增加或减少的数目绝对值),如果满足,那么,该簇将不再参与下一轮迭代调整。针对航运信息数据,通过多次Monte-Carlo仿真实验测试,在该算法中,通常设置阀值在00.1之间,阀值在00.05之间。阀值可以只设置任意一个,也可以两个都设置。但是,值得注意的是,如果阀值设置过大,会严重影响聚类精度,因此,要合理地设置阀值的大小。这样,有了调整阀值条件后,参与迭代的样本数目大大减少,从而可以节约大量的运算时间.算法的基本流程如下:输入:聚类个数k,数据集X,阀值和输出:k个簇集V和簇中心C：初始化簇中心;可以使用传统的k-means随机初始化法,也可以采用k-means初始化法。：设最终簇集和中心都为空集V=和C=,迭代次数t=1。计算每个对象与中心的距离,按最小距离原则将数据划分到相应的簇,形成k个簇,并记录各簇中样本个数,设置初始的调整标号.：重新计算每个簇中心和各个簇中心移动的初始位移,其中,表示簇中样本的个数.:计算标号的簇中每个样本对象与参与调整簇的中心的距离,按最小距离原则将数据划分到相应的簇,形成k个新簇,并记录各簇中样本个数,求得调整的样本个数.：重新计算参与调整簇中心和簇中心移动的位移.：如果簇的中心满足或者,那么,将簇中的调整标号设置为0,即=0,并将和添加到最终簇及簇中心集合中,即和.：如果形成了包含k个簇的集合,即V=k,那么,终止迭代,算法结束,否则,迭代次数t=t+1,跳转到步骤4.聚类结果的评价主要两个方面:聚类精度和算法时间复杂度.算法性能分析实验如下:(1)聚类精度分析传统的k-means算法,所求的解往往是局部最优解,使用k-means算法在一定程度上减少了局部最优可能性,但从聚类过程看,求的解也往往是局部最优解,本文提出的调整簇阀值的加速算法,由于采用了近似的k-means算法,导致该算法只能找到局部近似解。适当地设置阀值,也能得到和k-means算法相当的聚类精度.(2)时间复杂度分析传统k-means算法的时间复杂度为,其中,n为样本个数,k为簇数目,d为数据维数,t为迭代次数.5.2.3、相似度计算二元相似性度量算法分析 传统的度量是采用基于两幅图像对应像素位置的位与位匹配的方法。但是在被动篡改检测中,是没有可以用来比对的原始图像的。所以,在基于二元相似性度量的图像被动篡改检测中,我们考虑比较图像位平面的二元纹理统计特征. 首先,我们用表示第i像素和它的K连通像素(或K邻居像素,K=4或K=8)的顺序关系。在本文中,我们选择K=4,表示选择N,W,S,E(上!下!左!右)四个连通像素,i表示图像所有像素点,假定图像的大小是M和N。我们定义模板函数如式(6)所示： （19）这样，对于一个中间像素点i，它与其四连通像素的关系可以表示为式： （20）对于中心像素，我们定义它的协议变量(agreement)为下式： （21）其中，是克罗内克函数(Kronecker delta function)，其定义如式： （22）因此对于,有： （23） 很显然，函数表示中心像素点和其连通像素点之间的转换关系。对于中心像素点，它与其四连通像素的关系可以重新表示为式： 累积一协议量(accumulated agreements)的定义如下所示: ， ， 这四个变量是一幅二进制图像的一阶伴随变量.5.2.4、二叉树搜索排列二叉树搜索排序定义：二叉树的定义：空或有一个根，根有左子树、右子树；而左右子树本身又是二叉树。如图所示：B图3：LCDFE1、二叉树的性质：性质1：在二叉树的第i层上至多有2i-1个结点；性质2：高度为k的二叉树至多有2k-1个结点；性质3：二叉树的叶子结点数n0等于度为2的结点数n2+1；性质4：具有n个结点的完全二叉树高度为log2n+1；性质5：对一棵有n个结点的完全二叉树按照从第一层（根所在的层次到最后一层，并且每一层都按照从左到右的次序进行编号。根结点的编号为1，最后一个结点的编号为 n.1：对任何一个编号为i的结点而言，它的左儿子的编号为2i( 若2i= n) ,而右儿子的编号为2i+1(若 2i +1 = n).2：对任何一个编号为j的结点而言，它的父亲结点的的编号为j/2。根结点无父结点.2、二叉树的遍历:设 N 代表根节点,L代表左子树,R代表右子树。a. 前序（或先序）：如果二叉树为空，则操作为空：否则访问根结点；前序遍历左子树；前序遍历右子树。记为：NLRb. 中序：如果二叉树为空，则操作为空：否则中序遍历左子树；访问根结点；中序遍历右子树。记为：LNR. c. 后序 ：如果二叉树为空，则操作为空：否则后序遍历左子树；后序遍历右子树；访问根结点。记为：LRN.B图4：CLDFBEXW3、二叉树的作用:1查找当前第K大元素；2查询某元素在当前所有元素中的排名. （24） 其中，是克罗内克函数(Kronecker delta function)，其定义如式： （25） 因此对于,有： 很显然，函数表示中心像素点和其连通像素点之间的转换关系.对于中心像素点，它与其四连通像素的关系可以重新表示为式： 累积一协议量(accumulated agreements)的定义如式所示 ， ，这四个变量是一幅二进制图像的一阶伴随变量.5.2.5、构建模型的计算机仿真求解通过基于MATLAB编程的算法求解分析我们得到了如下表格的碎纸图片序号排列矩阵：表4: 附件三的碎纸片排列矩阵0123456104905406514318600220610190780670690993168100076062142030403814804616102403550711560831322000176014128003159082199709403408418309004781250131821091970169029064111201005092100072081381581260681108914616215411404078910111205719217811819009516209613107906311604102314719105017908118912210313019308003302219801513313501207316020316912104212414407711218411018706610615018004803707505504417504514700013705315120715514018510813141516171819011022129028091188141163072006177020052036120086195026001087018088167025008009105074170205085152165027060134039031051017115176149097136164127058043021170157181204139145206010104098172171059056093153070166032196117004101113194119123表5：对于附件四的碎纸片排列矩阵0123456119107501115419018422011481701961980943086051107029040158401919409314108812151591390011290632386020041108116136073720802100704906111980700840600140681749132181095069167163101710420662050101571108107712820013105278910111213002104180064106004149113164078103091080101186098024117150005059126105155114176182151153053038123120175085036207135015076043199033142168062169054192137195008047172156096166188111144206003130074145083134055018056125140193087089048072141516171819032204065039067147026100006017028146058092030037046127022057202071165082050160187097203031045173079161179143133118189162197112023099122090185109034013110025027178035016009183152044012177124000102115（附件3和4的碎纸片拼接结果详见附件）5.3、问题三模型的建立及求解5.3.1、蚁群算法实现图像最优化比对搜索 蚁群算法（ACA)是受到蚂蚁群体寻找食物行为的启发而提出的一种基于蚁群的模拟进化算法。一般来讲群体随机搜索算法常用于解决特定的组合优化问题。自然界中蚂蚁搜索食物过程是一个不断聚类的过程，食物就是聚类中心，其主要思想是:将每个数据看作一个蚂蚁，蚂蚁分别聚集到个聚类中心，到的距离为，采用欧式距离计算: （26）其中:m表示每个数据特征的维数;P加权因子，根据像素分量对信息测度影响程度设定。 假设算法中的蚂蚁具有一定的记忆能力，能根据2幅待拼接图像上互信息浓度选择转移方向，从而引导蚂蚁向互信息测度最大值的方向移动. 蚁群算法处理的问题一般具有如下特点:搜索空间是离散的;有一组有限的约束条件;有一个代价函数，为搜索算法生成的解计算对应的代价，解的每一部分都会对解的代价产生影响;有一个有限的节点集合，用于构建解;有一个有限的节点间的可能转移的集合;有一个节点序列的有限集合，用于表示所有的有效组合，来定义完整的搜索空间，组合的有效性、可行性由约束条件决定.可以将图像拼接过程中图像A中的窗口在B搜索迭代的过程看成是“只蚂蚁的城市模型”，从而使得蚂蚁可以在这个模型图上进行爬行，并且保证了蚂蚁爬行所获得的每个哈密顿路都对应一个互信量最大的组合. 蚁群算法的参数进行如下设置: 为信息素的相对重要程度; 为启发式因子的相对重要程度; 为信息素蒸发系数，表示信息素的持久性系数，且 （27）其中:为窗口信息素含量;为启发函数;N为循环遍历次数.5.3.2、基于ACA搜索的碎纸片图像SIFT特征点匹配基于图像特征的方法，首先要对待配准的两幅图像进行处理，提取满足特定应用要求的特征集，然后将这两组特征集进行匹配对应，生成一组对应特征对集，最后利用这组特征对之间的对应关系估计出全局变换参数。基于图像特征的方法，在特征提取后得到的特征点的数量将会大大减少，因此可以提高配准的速度，但其配准的效果很大程度上还取决于特征点的提取精度以及特征点匹配的准确度。基于图像特征的配准方法主要困难在于如何提取和选择鲁棒的特征，以及如何对特征进行匹配，其中要克服，由于图像噪声和场景中出现遮挡现象所引起的误匹配的问题。常用的图像匹配特征有点、直线、曲线等.目前大多数文献都是采用点特征进行图像之间的配准，从而实现图像的拼接。基于图像特征配准的方法的主要优点是它提取了图像的显著特征，大大压缩了图像的信息量，计算量小，速度较快，而且它对图像灰度的变化具有鲁棒性。但另一方面，正是由于只有一小部分的图像灰度信息被使用了，所以这种方法对特征提取和特征匹配的错误更敏感，需要可靠的特征提取和特征一致性.图5：基于特征点的图像拼接流程图 该问题总结的算法的搜索过程可以理解为:在2幅或多幅图上利用蚁群作为搜索窗口，通过互信息测度的方式代替搜索出的2幅图像之间在目标、背景、边界和噪声等内容中的特征向量的过程。依据互信息量最大作为搜索并且不断改变的方向，最终朝着使互信息量最大的方向搜索，从而确定出2幅图像之间存在的匹配区域.碎片全局拼接如果有两个碎片证明是匹配的，那么把匹配之后的碎片组成一个新的人工碎片，重复以上的步骤，直到所有的碎片拼接完毕。用这种方式，拼接的过程形成了岛屿状。计算每个碎片的轮廓长度。假如M是旋转链的轮廓像素的总数，N是固定链的轮廓像素的总数。那么旋转链每次从#（M-j），j=1,2,.,M-1绕着固定链# (N-L),L=1,2,. . . .N-1进行比较，直到完成匹配.5.3.3、模型三的MATLAB求解表6：对于附件五的碎纸片排列矩阵第一面01234561078b111b125a140a155a150a2089a010b036a076b178a044a3186b153a084b042b030a038a4199b011b161b169b194b173b5088b107a149b180a037b191a6114a184b179b116b207a058a7146a171b031a201a050a190b8165b195a118a157a168a046a9003b007b085b148b077a004a10023b133a048a051b095a160b11099a043a096b109a123a006a78910111213183b174b110a066a108a018b029a025b192a124b022a120b144a079a121a098a094b061a137b045a138a206b156a034a181b198b087a132b065b115b166b001b151b170b041a158a197a154b028b012a017b120b092b019b016b177b053b202a021b067a063b075b167a117b008b068b069a032a074b126b176a158a000b119a033b071b052a062a129b118b104a134a113a026b049b091b106b141516171819198b081b164b020a047a136b014a059a060b147a152a005a056b131b187b086b200b143b093a072b175a097a039b083a070b139b002a162b203b090a064b208a142a057a024a013a130a163a193b073b159a035a188a127a040a182b122a172a180b027a135b141a204b105a101a015b205a082b145a009b100b055b103a112a196b054b表7：对于附件五的碎纸片排列矩阵第二面01234561136a047b020b164a081a198a2005b152b147b060a059b014b3143a200a086a187a131a056a4083b039a097b175b072a093b5090b203a162a002b139a070a6013b024b057b142b208b064a7035b159b073a193a163b130b8172b122b182a040b127b188b9105b204a141b135a027b180a10009a145b082a205b015a101b11054a196a112b103b055a100a78910111213029b018a108b066b110b174a183a079b144b120a022b124a192b025a138b045b137a061b094a098b121b132