岩石颗粒图像中的斜长石识别

周磊，滕奇志，何海波

（1.四川大学电子信息学院图像信息研究所，成都 610065；2.成都西图科技有限公司，成都 610024）

0 引言

在石油地质领域中，研究人员通过偏光显微镜观察岩石的光学特性和表面纹理特征，从而对不同岩石进行鉴别和成分分析，但在实际应用中岩石样本目标颗粒众多，人工鉴定工作量大，因此利用计算机实现岩石矿物成分分类识别具有重要研究价值。针对斜长石类矿物易呈现出相互平行的特殊纹理直线的特征，通过提取这些具有一定规律的直线特征达到对斜长石族矿物及其他岩石颗粒进行区分的目的。

直线特征的检测与提取一直以来都是图像处理和分析的一个重要问题，其中经典的算法包括Hough变换和LSD(Line Segment Detector)局部直线检测算法。Hough[1]变换是传统的直线检测方法之一，具有较强的抗噪声能力并能够准确给出所提取直线段特征参数。然而，将Hough变换应用于岩石薄片图像时，在密集点处出现大量共线点，反映在Hough的参数空间中具有许多伪峰点，而在实际的应用中伪峰点在图像上无法找到其对应的直线，从而造成误提取。

LSD局部直线检测算法是由VGR Grompone等[3]提出，算法首先计算每个像素点的梯度信息，并根据区域增长生成直线段候选区域，再通过外接矩形计算出区域所代表的直线段，并选择其长轴作为候选直线段区域。该算法的特点是提取速度快，无需调整任何参数。不少学者在经典LSD算法基础上进行了研究，参考文献[4]提出了一种基于区域生长，通过图像梯度信息的局部直线段检测方法，参考文献[7]对参考文献[4]做了改进，减少了计算量，提高了检测率，但是该算法采用了固定阈值大小，针对不同应用需要相应的调整其大小。文献[6]基于LSD局部直线检测算法，提出直线区域有效性检测策略的直线段检测算法，文献[11]参考文献[6]使用相同的直线区域有效性检测策略的直线段检测算法，提出了一种基于动态阈值的区域生长的LSD算法，可以降低方向准确度使得候选区域更加连续，从而获得更准确的直线。

图3 梯度模板

我们将LSD直线检测算法应用在颗粒目标上时发现，由于局部检测算法自增长的特点，颗粒中长线段局部模糊或被遮挡的原因，经常被割裂为多条直线，造成大量噪声和断裂现象，该现象的出现往往会造成斜长石颗粒中的直线特征的误识别。针对此问题，本文提出了一种基于LSD候选直线段分段拟合的方法。该算法既能够实现精准、快速的直线段特征提取，又能提高复杂纹理目标下检测直线的有效性和直线段的完整性。

1 岩石薄片纹理特征选择和提取

1.1 岩石薄片纹理特征的选择

图像特征的提取是图像识别和分类的前提，根据图像的特点通过某种算法有效地区分出图像中不同类型的物体的特征向量或者特征参数。一般来说，特征提取分为基于形状特征提取、空间特征提取、颜色特征提取和纹理特征提取4种类型。每种类型都有不同的方法，针对不同的图像，每种方法呈现出较大差异，故应结合实际图像来选取正确的特征提取方法。本文主要针对岩石矿物薄片在偏光镜下易呈现的双晶、解理、裂理等纹理的直线特征，通过提取这些具有一定规律的纹理特征可以有效的完成对斜长石与其他岩石矿物的区分。

下图1分别介绍了三种岩石矿物在偏光显微镜下的纹理特征，（a）石英矿物图像，其内部特征只有少量的，离散的噪声，纹理细节较少，颗粒表面整洁光滑；（b）斜长石矿物图像，其特征表现为由于双晶、糙面等现象使得颗粒的纹理细节比较清楚，纹理比较丰富，且存在若干条一定间距的平行纹理直线；（c）碱性长石矿物图像，同样存在大量的纹理信息，但该纹理信息不具有规律性，既不满足斜长石的纹理特征，同时跟石英相比具有大量的无规则信息。

图1 石英、斜长石以及碱性长石矿物纹理特征的区别

1.2 基于LSD候选直线段分段拟合算法的特征提取

复杂岩石颗粒图像中的斜长石识别算法主要分为图像预处理、基于LSD直线段的候选区域提取，之后根据提供的候选直选段对其进行二维特征的相似度进行聚类，最后根据图像直线密度进行有效性的检测并进行斜长石识别。对岩石颗粒图像整体识别过程流程图如图2所示。

（1）图像的梯度

图像的梯度反映了图像中相邻像素间的灰度值的变化关系，计算图像的梯度对于检测图像的边缘结构有着十分重要的意义和作用。对于数字图像I中像素点（x,y），i（x,y）为该像素的灰度值，gx、gy是图像坐标系的两个梯度方向，对于不同的梯度算子，模板的大小也会直接影响图像的细节信息，模板过大容易导致图像细节丢失，模板过小虽然节省了计算量，但对噪声敏感度增加，增加不必要的冗余信息。本文参照文献[4]采用最简单的 2×2模板，通过式（1）-（4）可计算出图像的梯度大小和方向。

梯度的模值：

梯度方向的计算：

梯度的模板：

（2）图像区域生长策略

针对得到的梯度图像信息，采用区域生长策略来获得具有相似性的连通区域。区域生长是一个不断迭代的过程，通过种子点的选择，根据一定的生长准则，把相邻像素点的相似属性作为分类的标准，逐步迭代分类，直至所有点都迭代完毕或者满足终止条件。因此可以通过该方法将具有相近梯度信息的像素点连通起来构成候选区域，该算法思想简单并且能够提供较好的边界信息。

本文参考文献[6]的贪心算法的区域生长策略，该策略既减少算法的时间代价，又提高线特征提取速度，同时经过生长区域得到的像素点只能属于一个区域，所以在生长的过程中需要对其标注不同区域的标签，直至所有的点遍历完毕。在数字图像中，梯度值较大的像素点对应着更加明显的边界，直线段存在边缘信息上的可能性也增加了。因此需对像素梯度的大小进行排序，本文采用文献[7]的伪序法，该方法是以线性时间为代价的桶排序思想的方法，使排序速度得到大幅提升。排序后将像素梯度值最大的点作为种子点进行基于梯度信息的区域生长，根据生长准则（邻域点的梯度方向与区域主体的梯度方向的差值在一定的阈值θ范围之内）标记出属于该区域的所有点，再将剩余的梯度的模值较大的点作为种子点，依次进行迭代，直至某像素点不符合生长准则或者遍历完所有像素点。如果阈值θ太小，会导致检测到的线段断裂，甚至部分直线段无法检测，如果阈值θ设置太大，区域生长的过程中会出现许多误连接，导致生长出的区域内包含太多噪声点，容易将噪声点误检测为直线段的一部分。

（3）基于候选区域的相似直线拟合

通过以上步骤得到直线段的候选区域后需要将该区域拟合成一条完整连续的直线段，本文采取了通过外接最小矩形来描述候选区域，计算得到的每个像素作为加权值，并参与外接矩形中心位置的计算，具体公式如式（5）、（6）所示：

其中R代表某个区域，G（j）代表R内第j点的梯度大小，x（j）是 R 内第 j点的 x坐标值，y（j）是 R 内第 j点的y坐标值。在计算矩形中心后，将矩形的主方向设置为直线主方向，将矩阵M中最小特征值的特征向量的方向作为矩形的方向，矩阵M如式（7）所示。

最后，在矩形能够包含区域内所有点的基础上，将矩形的长轴作为候选目标直线段。

由于LSD算法本身的局限性，目标颗粒中同一条解理线的直线段之间存在严重的断裂现象，导致在识别判断的过程对同一条解理线进行多次重复计算，从而造成整体目标的误识别。为了获得较完整直线段数据，本文利用直线段的距离特征对具有相似倾斜角的线簇进行相似拟合，将倾斜角相近且距离相近的直线簇拟合为一条直线，这样做的主要目的是确保复杂纹理图像中直线段的完整性，为后期的统计分类提供可靠的数据信息，提高斜长石的分类的准确率。其主要步骤如下。

步骤1：首先对直线段进行一定的预处理操作，删除冗余直线段。在复杂纹理薄片直线检测的过程中，针对线段长度过短或者总体数量较少的情况，将这些直线段视为噪声线段进行过滤，其主要目的是为了获得薄片图片中更加明显的长直线特征。

步骤2：笛卡尔坐标转极坐标。将直线倾斜角θ和直线到原点的距离d，用极坐标表示。

根据公式：

斜截式：截距为b,斜率为k,倾斜角为θ

点到直线的距离公式:

根据公式（8）-（10）求得描述直线的二维特征：倾斜角θ和直线到原点（图片的左下角作为原点）的距离d。

如图4所示，图a为碱性长石的样本颗粒，图b是基于LSD算法下检测出来的直线纹理特征，图c表示统计出的每条直线的倾斜角度信息，从该直方图可以发现直线纹理的倾斜角度主要集中在125°左右，图d表示的是每条直线到原点的不同距离，根据纹理直线的倾斜角度和直线到原点之间的距离的直线二维特征，寻找倾斜角度接近和直线到原点距离相近的直线并归类。

图4 统计去噪后的直线二维特征直方图

步骤3：直线段的描述。在极坐标下，倾斜角θ和原点距离d可以描述不同的直线段，每条直线段对应唯一的倾斜角θ和原点距离d，对该二维特征分别进行直方图统计。已知斜长石的特征主要是存在一定间距的平行直线：对于平行直线，倾斜角相近，为了保证平行直线在的倾斜角误差范围内，在数字图像处理技术中，若有两像素点的梯度方向之差direction（g（x1，y1）。g（x2，y1））<2π/n，则可认为在方向准确度为 n 时两像素的梯度方向一致。根据神经感知学理论[8]以及数值实验，通常n取8～16，此处根据实验效果设置n为72，其步长角度为5°，因为步长角度过大，会造成角度过大的分为一簇，若步长角度太小，直线簇类过多，都会影响实验精度。

步骤4：分类拟合。在已知相互平行的直线簇中，根据距离阈值D，依次选取其中一条并与其他直线距离比较获得对应的差值距离Δd，若Δd=(d1-d2）D的情况，循环结束。将拟合后的新直线段，代替被拟合的多条直线段，如下图5所示，图b描述的是需要待拟合的直线段，图c是将拟合后的直线段整合到剩余直线段中。

图5 多条相似直线段的拟合过程

步骤5：直线密度有效性判别。k为图像中主方向直线的条数。若k>D，则图像张检测的直线密度判定有效，可作为斜长石判定的候选图像，否则，该岩石颗粒判定为非斜长石。本文设置D=2，直线密度有效性判别将会解决因为存在较少的平行直线单元情况而误判别为斜长石的状况。

（4）斜长石的识别

经过拟合后的直线，斜长石颗粒直线的纹理特征更加明显。针对上述过程中获得的直线倾斜角θ和原点距离d，考虑其中多条直线平行度，便可判定该岩石颗粒是否属于斜长石颗粒。

根据公式：

标准特征近似参数γ：num为主方向直线的条数，N为所检测直线的总数：

根据公式（8）-（11），若存在γ>C1（0

2 实验结果及分析

本文算法在检测复杂岩石纹理颗粒直线特征时，其主要目的是从各类岩石颗粒中识别出斜长石的纹理特征，因此本文算法针对该问题，提出了两项要求：1）对直线的检测过程必须要保证直线特征的直观性和有效性，即得到高质量的斜长石的纹理特征；2）滤除冗余噪声直线，在复杂的岩石薄片图像中，夹杂中许多的噪声信息，为了保证纹理特征的有效提取，尽可能减少因为噪声直线，对提升斜长石的识别率有很大的帮助。

从图7可以观察到两种基于LSD直线算法，图a和e是斜长石样本颗粒，图b和f是传统LSD直线检测算法，存在大量噪声直线段以及直线断裂现象，因为图中存在非平行直线段和噪声直线的干扰，在此基础上用相同的斜长石判别方法进行判断，其结果属于非斜长石类，图c和g是文献[11]算法，该算法基于LSD算法，采用了动态阈值的区域生长，能够更好地抑制纹理边缘所造成的伪直线，减少了直线段中的断裂现象,图d和h是本文算法，从图中可以观察到相比图b、f、c和g，本文算法滤除了大量冗余直线段，同时成功拟合了部分断裂直线段，呈现出具有一定规律的平行直线段并且线条完整与实际纹理特征相符合，达到提供更加完整连续平滑直线段的目的，并且克服了LSD局部直线检测中出现的断裂直线段缺点的问题，为斜长石直线特征的识别提供了有力的依据。

为了验证该算法的识别率，本文采用了三组不同的岩石样本实验数据，第一组为斜长石颗粒样本，第二组为石英颗粒样本，第三组是碱性长石和石英的混合颗粒样本，对三组不同的测试样本检测出斜长石颗粒的正确识别率以及误识别率，测试结果如下表所示。表1可以明显观察到采用本文算法可以对岩石颗粒中的斜长石具有较好的识别能力，提升了正确识别率的结果，同时从表2可以观察到本文算法对碱性长石和石英有很好的区分能力，降低了误识别率的结果。

表1 本文算法针对斜长石样本颗粒识别结果

图7 文献[11]直线检测算法以及本文算法对比效果图

表2 本文算法针对石英和碱性长石样本颗粒的检测结果

3 结语

本文提出了一种基于LSD直线候选区域根据相似性拟合的算法。该算法主要通过图像的梯度信息进行区域生长，获得具有相似梯度区域的连通区间，然后用外接最小矩形的长轴描述候选直线段，再对候选直线段进行相似拟合，获得能够表征多段分裂直线段的直线段结构。通过实验和理论分析表明，本算法具有以下主要特点：1）算法具有很强的鲁棒性，能够在复杂的纹理图片中检测到直线；2）算法检测的直线段更加精确。相比Hough变换中存在的大量误检和LSD直线检测过多断裂直线段，该算法保证了获得复杂纹理图片中连续直线的完整性，有效地抵抗了细小的边缘结构区域；3）该算法计算量小，继承了文献[7]采取的伪排序算法滤除了复杂纹理中的大量噪声，节省了计算量。

[1]村上，伸一.Use of the Hough Transform to Detect Lines and Curves in Pictures[C].Commun.ACM 1972.1972:11-15.

[2]Galamhos C，Matas J，Kittler J.Progressive Probablistic Hough Transform for Line Detection.Proceedings of Computer Vision and Pattern Recognition，1999：23-25.

[3]Kultanen P,Xu L,Oja E.Randomized Hough Transform[C].Proc.Of IEEE Pattern Recognition,1990,631-635

[4]Rafael Grompone von Gioi,Jérémie Jakubowicz,Jean-Michel Morel,and Gregory Randall,LSD:a Line Segment Detector,Image Processing on Line,2（2012）：35-55.

[5]夏军营,徐小泉,熊九龙.利用梯度信息快速提取直线边缘特征[J].中国图象图形学报,2012,17（8）:987-994.

[6]董银文,苑秉成,石钊铭,等.基于直线特征的航拍图像机场跑道自动识别算法[J].系统工程与电子技术,2013,35（4）:876-879.

[7]覃勋辉,马戎,付维平,等.一种基于梯度的直线段检测算法[J].光子学报,2012,41（2）.

[8]Desolneux A,Moisan L,Morel J M.Computational Gestalts and Perception Thresholds[J].Journal of Physiology-Paris,2003,97（2-3）:311-324.

[9]Xu Lei，Oja E.Randomized Hough Transform（RHT）:Basic Mechanism,Algorithms,and Computational Complexities[J].CVGIP:Image Understanding,1993,57（2）:131-154.

[10]Lee D,Park Y.Discrete Hough Transform Using Line Segment Representation for Line Detection[J].Opt Eng,2011,50（8）:087004.

[11]闵涛,蒋宜勤,滕奇志,等.矿物显微图像中的直线段特征提取算法[J].微型机与应用,2016,35（6）:35-38.