CN110827306A - 基于自适应密度峰值聚类的紫色土图像分割提取方法 - Google Patents

基于自适应密度峰值聚类的紫色土图像分割提取方法 Download PDF

Info

Publication number
CN110827306A
CN110827306A CN201910991121.XA CN201910991121A CN110827306A CN 110827306 A CN110827306 A CN 110827306A CN 201910991121 A CN201910991121 A CN 201910991121A CN 110827306 A CN110827306 A CN 110827306A
Authority
CN
China
Prior art keywords
matrix
point
imgin
image
imggrays
Prior art date
Legal status (The legal status is an assumption and is not a legal conclusion. Google has not performed a legal analysis and makes no representation as to the accuracy of the status listed.)
Granted
Application number
CN201910991121.XA
Other languages
English (en)
Other versions
CN110827306B (zh
Inventor
曾绍华
唐文密
詹林庆
罗达璐
Current Assignee (The listed assignees may be inaccurate. Google has not performed a legal analysis and makes no representation or warranty as to the accuracy of the list.)
Chongqing Normal University
Original Assignee
Chongqing Normal University
Priority date (The priority date is an assumption and is not a legal conclusion. Google has not performed a legal analysis and makes no representation as to the accuracy of the date listed.)
Filing date
Publication date
Application filed by Chongqing Normal University filed Critical Chongqing Normal University
Priority to CN201910991121.XA priority Critical patent/CN110827306B/zh
Publication of CN110827306A publication Critical patent/CN110827306A/zh
Application granted granted Critical
Publication of CN110827306B publication Critical patent/CN110827306B/zh
Active legal-status Critical Current
Anticipated expiration legal-status Critical

Links

Images

Classifications

    • G—PHYSICS
    • G06—COMPUTING; CALCULATING OR COUNTING
    • G06T—IMAGE DATA PROCESSING OR GENERATION, IN GENERAL
    • G06T7/00—Image analysis
    • G06T7/10—Segmentation; Edge detection
    • G06T7/13—Edge detection
    • G—PHYSICS
    • G06—COMPUTING; CALCULATING OR COUNTING
    • G06F—ELECTRIC DIGITAL DATA PROCESSING
    • G06F18/00—Pattern recognition
    • G06F18/20—Analysing
    • G06F18/23—Clustering techniques
    • G—PHYSICS
    • G06—COMPUTING; CALCULATING OR COUNTING
    • G06T—IMAGE DATA PROCESSING OR GENERATION, IN GENERAL
    • G06T7/00—Image analysis
    • G06T7/10—Segmentation; Edge detection
    • G06T7/11—Region-based segmentation
    • G—PHYSICS
    • G06—COMPUTING; CALCULATING OR COUNTING
    • G06T—IMAGE DATA PROCESSING OR GENERATION, IN GENERAL
    • G06T7/00—Image analysis
    • G06T7/90—Determination of colour characteristics
    • G—PHYSICS
    • G06—COMPUTING; CALCULATING OR COUNTING
    • G06T—IMAGE DATA PROCESSING OR GENERATION, IN GENERAL
    • G06T2207/00—Indexing scheme for image analysis or image enhancement
    • G06T2207/20—Special algorithmic details
    • G06T2207/20004—Adaptive image processing
    • G—PHYSICS
    • G06—COMPUTING; CALCULATING OR COUNTING
    • G06T—IMAGE DATA PROCESSING OR GENERATION, IN GENERAL
    • G06T2207/00—Indexing scheme for image analysis or image enhancement
    • G06T2207/30—Subject of image; Context of image processing
    • G06T2207/30181—Earth observation

Landscapes

  • Engineering & Computer Science (AREA)
  • Theoretical Computer Science (AREA)
  • Physics & Mathematics (AREA)
  • Computer Vision & Pattern Recognition (AREA)
  • General Physics & Mathematics (AREA)
  • Data Mining & Analysis (AREA)
  • Artificial Intelligence (AREA)
  • Bioinformatics & Cheminformatics (AREA)
  • Bioinformatics & Computational Biology (AREA)
  • Evolutionary Biology (AREA)
  • Evolutionary Computation (AREA)
  • Life Sciences & Earth Sciences (AREA)
  • General Engineering & Computer Science (AREA)
  • Image Analysis (AREA)

Abstract

本发明提供一种基于自适应密度峰值聚类的紫色土图像分割提取方法，所述方法包括以下步骤：S1：对含有紫色土区域的紫色土彩色图像进行可分离性灰度变换，获得灰度图像Ⅰ；S2：利用自适应密度峰值聚类算法对灰度图像Ⅰ进行初步分割，获得初步分割后的二值图像Ⅱ；S3：对二值图像Ⅱ进行边界提取处理，获得紫色土壤区域的边界矩阵；S4：对提取的边界矩阵进行填充，获得二值图像Ⅲ；S5：求出二值图像Ⅲ与所述含有紫色土区域的彩色图像的哈达玛积，得到只包含有紫色土图像的分割图像。本申请可将紫色土区域图像从背景中准确、完整分割出来，且在分割过程实现紫色土的自适应的分割，具有分割速度快、准确、完整的有益技术效果。

Description

基于自适应密度峰值聚类的紫色土图像分割提取方法
技术领域
本发明涉及图像分割提取方法，尤其涉及一种基于自适应密度峰值聚类的紫色土图像分割提取方法。
背景技术
随着人工智能技术的不断发展，其已广泛应用于生产、生活的诸多领域；在人工智能领域的众多分支中，机器视觉得到快速发展，简单说来，机器视觉就是用机器代替人眼来做测量和判断。在农业生产中，土壤的土种及其土属的识别完全依赖于农业专家，众所周知，土壤的识别在农业生产中是非常重要的，因土壤识别对识别人的专业技能要求高，如现有土壤的土种及其土属的识别完全依赖于农业专家，阻碍了土壤识别的普及。但随着现代农业的发展，以及人工智能技术的不断发展，使机器视觉来识别土壤成为可能。机器视觉辨识土壤是对野外自然条件下拍摄的具有复杂背景的土壤图像进行识别并精准判断土壤所属土种，但野外自然环境下拍摄的紫色土彩色图像一般都包含有农作物、地衣、苔藓、杂草等复杂背景，会对机器视觉识别造成干扰。如何排除背景区域的干扰，精准地从含有土壤区域的彩色图像中提取出土壤区域图像来是机器视觉识别紫色土土种的基础工作，同时也是目前亟待解决的技术难点。紫色土是我国西南地区最主要耕地之一，紫色土是机器视觉土壤识别的研究重点。
因此，亟需一种能快速、准确地将紫色土土壤区域图像从含有紫色土区域的彩色图像中分割提取的方法。
发明内容
有鉴于此，本发明提供一种基于自适应密度峰值聚类的紫色土图像分割提取方法。
本发明提供一种基于自适应密度峰值聚类的紫色土图像分割提取方法，其特征在于：包括以下步骤：
S1：对含有紫色土区域的紫色土彩色图像进行可分离性灰度变换，获得灰度图像Ⅰ；
S2：利用自适应密度峰值聚类算法对灰度图像Ⅰ进行初步分割，获得初步分割后的二值图像Ⅱ；
S3：对二值图像Ⅱ进行边界提取处理，获得紫色土壤区域的边界矩阵；
S4：对提取的边界矩阵进行填充，获得二值图像Ⅲ；
S5：求出二值图像Ⅲ与所述含有紫色土区域的彩色图像的哈达玛积，得到只包含有紫色土图像的分割图像。
进一步，所述步骤S1包括：
S101：将表征含有紫色土区域的紫色土彩色图像的矩阵ImgIn从RGB颜色空间转换到Lab颜色空间，获得矩阵ImgIn'，基于Lab颜色空间计算图像相似度矩阵Sx,y，以及矩阵ImgIn'中点ImgIn＇x高相似度点集合和低相似度点集合；
S102：在Lab颜色空间中建立灰度变换优化模型，采用离散搜索求解该优化模型，获得灰度矩阵ImgGray，即灰度图像Ⅰ。
进一步，所述步骤S101包括：
S1011：在Lab颜色空间中，获取L、a、b分量的离散二维信息熵，所述L、 a、b分量的离散二维信息熵联立式子(1)至(3)确定，
Figure RE-GDA0002360084660000021
Figure RE-GDA0002360084660000022
c＝{L,a,b} (3)
其中，Hc表示单分量图像的离散二维熵，
Figure RE-GDA0002360084660000033
表示中间变量，
Figure RE-GDA0002360084660000032
表示单分量图像像素灰度值d与像素8邻域灰度均值e构成的二元组(d,e)出现的频数， eight表示图像的高度，width表示图像的宽度；
S1012：计算矩阵ImgIn'的相似度矩阵Sx,y，所述相似度矩阵Sx,y采用如下方法确定：
Figure RE-GDA0002360084660000031
其中，Sx,y表示矩阵Img'I的相似度矩阵，x＝1,...,width*height， y＝1,...,width*height，Hc表示单分量图像的离散二维熵，ImgIn＇c,x表示矩阵ImgIn'中点ImgIn＇x的c分量的值，ImgIn＇c,y表示矩阵ImgIn'中点ImgIn＇y的c分量的值；
S1013：计算矩阵ImgIn'中点ImgIn＇x高相似度点集合和低相似度点集合，
所述高相似度点集合采用如下方法确定：
set_highx＝{ImgIn＇y|Sx,y≥meanSx} (5)
其中，set_highx表示点ImgIn＇x的高相似度点集合，ImgIn＇y表示矩阵ImgIn'中除ImgIn＇x外的其他点，Sx,y表示矩阵ImgIn'的相似度矩阵，meanSx表示点ImgIn＇x与所有点间相似度的均值，
所述低相似度点集合采用如下方法确定：
set_lowx＝{ImgIn＇y|Sx,y＜meanSx} (6)
其中，set-lowx表示点ImgIn＇x的低相似度点集合，ImgIn＇y表示矩阵ImgIn'中除ImgIn＇x外的其他点，Sx,y表示矩阵ImgIn'的相似度矩阵，meanSx表示点ImgIn＇x与所有点间相似度的均值。
进一步，所述步骤S102包括：
S1021：建立优化模型计算Lab颜色空间亮度分量L的权重wL，颜色分量 a的权重wa,颜色分量b的权重wb,所述优化模型如下：
Figure RE-GDA0002360084660000041
其中，wL表示亮度分量L的权重，wa表示颜色分量a的权重，wb表示颜色分量b的权重，mheight表示图像的高度，nwidth表示图像的宽度,λ表示调节参数λ∈(0，G_hx和G_lx表示中间变量，
Figure RE-GDA0002360084660000043
其中，ImgIn＇y表示矩阵ImgIn'中除ImgIn＇x外的其他点，set_highx表示点ImgIn＇x的高相似度点集合，set_lowx表示点ImgIn＇x的低相似度点集合，Sx,y表示矩阵ImgIn'的相似度矩阵， ImgGrayx表示任意点ImgIn＇x灰度变换后的值，ImgGrayy表示任意点ImgIn＇y灰度变换后的值，L表示Lab颜色空间中的亮度分量，a表示Lab颜色空间中的颜色分量，b表示Lab颜色空间中的亮度分量，wc表示单分量的权重，ImgIn'c,x表示点 ImgIn＇x的某单分量值；
S1022：离散搜索求解优化模型(7)，获得紫色土灰度矩阵ImgGray，对所述紫色土灰度矩阵ImgGray采用如下方法获得：
a：求解临时矩阵ImgGrayTemp，
ImgGrayTempx＝wL*ImgIn＇L,x+wa*ImgIn＇a,x+wb*ImgIn＇b,x (8)
其中，ImgGrayTempx表示临时矩阵ImgGrayTemp中点x的值，wL表示亮度分量L的权重，wa表示颜色分量a的权重，wb表示颜色分量b的权重，ImgIn'L,x表示点ImgIn＇x的L分量值，ImgIn'a,x表示点ImgIn＇x的a分量值，ImgIn'b,x表示点ImgIn＇x的b分量值；
b：对矩阵ImgGrayTemp中元素向下取整，获得紫色土灰度矩阵ImgGray,即为灰度图像I。
进一步，所述步骤S2包括：
S201：计算点ImgGraySi的密度ρi和点ImgGraySi的距离δi，其中，ImgGraySi表示数据集ImgGrayS中第i个点；
S202：自助确定聚类中心；
S203：利用自适应密度峰值聚类算法分割紫色土灰度矩阵ImgGray，获得二值紫色土初分割矩阵ImgIniSeg，即二值图像Ⅱ。
进一步，所述步骤S201包括：
S2011：统计紫色土灰度矩阵ImgGray的一维直方图，将ImgGray的灰度级作为待聚类数据集ImgGrayS＝{ImgGrayS1,…,ImgGraySi…,ImgGraySN},包含N个数据， N＝max(ImgGray)-min(ImgGray),min(ImgGray)表示频数不为0的最小灰度级， max(ImgGray)表示频数不为0的最大灰度级；
S2012：计算点ImgGraySi的密度ρi，并记录ρi序排序时对应的原序号orderρi，其中密度ρi采用如下方法确定，
Figure RE-GDA0002360084660000051
其中，ρi表示点ImgGraySi的密度，N表示聚类数据集ImgGrayS的元素个数，frequencyj表示ImgGraySj的频数，dij为ImgGraySi与ImgGraySj间的欧式距离；
S2013：计算点ImgGraySi的距离δi，记录比ImgGraySi密度大的最近相邻点 ni，并将ImgGraySj赋值给ni，即ni＝ImgGraySj，对于密度最大的点，其ni＝0，其中，距离δi采用如下方法确定：
其中，δi表示ImgGraySi的距离，dij为ImgGraySi与ImgGraySj间的欧式距离，δj表示非点ImgGraySi的点ImgGraySj的距离。
进一步，所述步骤S202包括：
S2021：对ρi和δi分别进行min-max归一化处理，获得ρi'和δi'；定义任意点ImgGraySi到以ρi'为横轴，以δi'为纵轴构成的二维矩阵的左下角原点的线段作为三角形的斜边Hi，定义ImgGraySi到矩阵副对角线的垂线段作为三角形的对边 OSi，定义垂足点到原点的线段作为三角形的邻边ASi，
Figure RE-GDA0002360084660000062
OSi＝sinβi×Hi (12)
ASi＝cosβi×Hi (13)
其中，δi'表示进行了min-max归一化处理ImgGrayiS的距离，ρi'表示进行了 min-max归一化处理点ImgGraySi的密度，Hi表示任意点ImgGraySi到以ρi'为横轴，以δi'为纵轴构成的二维矩阵的左下角原点的线段作为三角形的斜边，OSi表示ImgGraySi到矩阵副对角线的垂线段作为三角形的对边，ASi表示垂足点到原点的线段作为三角形的邻边，βi表示Hi与副对角线的夹角βi＝|αi-45|，Hi与横轴的夹角αi＝arctan(δi/ρi)/π*180；
S2022：计算中心决策度量CDMi，其中，中心决策度量CDMi采用如下方法确定，
CDMi＝Hi-OSi (14)
其中，Hi表示三角形的斜边，OSi表示三角形的对边；
S2023：对CDMi采用最大类间方差法确定分类阈值T，当CDMi≥T时，对应灰度级点即聚类中心点；
S2024：初始化每个点的类标记Ci，若ImgGraySi是聚类中心点，ci＝k， k表示第k个类，否则ci＝-1；
S2025：按orderρi的顺序将非聚类中心点归类到比其高密度点的最近邻点 ni所属的类中：若ci＝-1，则ci＝cni，Ci表示数据集ImgGrayS中第i个点 ImgGraySi的类标记，cni表示比ImgGraySi密度大的最近相邻点，其计算见步骤 S3013；
S2026：创建一个与矩阵ImgGray大小相同的零矩阵ImgIniSeg，遍历ImgGray全部元素，若ImgGraym,n＝ImgGraysi，则ImgIniSegm,n＝ci，m为行号，n为列号；
S2027：在矩阵ImgGray中心点附近随机抽取7个3×3大小子阵，剔除均值最大和最小的子阵，计算剩余子阵全部元素的均值μ作为土壤中心土块均值；计算每一类的聚类中心与土壤中心土块均值μ的距离，距离最小的那一类，将矩阵ImgIniSeg对应位置元素置1，否则置0，获得初分隔矩阵ImgIniSeg，即为初分割二值图像II。
进一步，所述步骤S3中的边界提取包括以下步骤：
S311：对初分割矩阵ImgIniSeg中，搜索边界起点：扫描初分隔矩阵ImgIniSeg，寻找边界起点ImgIniSeg(m,n)，所述边界起点ImgIniSeg(m,n)需同时满足： ImgIniSeg(m,n)＝1且在ImgIniSeg(m,n)的4-领域内，ImgIniSeg(m,n)的右方或下方的点至少有一个的值为1；
S312：搜索整条边界并记录该条边界最小行号mmin与最大行号mmax，以及该条边界最小列号nmin和最大列号nmax；
S313：若最小行号mmin与最大行号mmax，以及最小列号nmin和最大列号nmax满足以下式子，若满足，则结束；若不满足，将矩阵ImgIniSeg和ImgBoundary 中行号在mmin和mmax之间且列号在号nmin和nmax之间值为1的点置为0，并进入步骤S311；
其中，mmin表示最小行号，mmax表示最大行号，nmin表示最小列号，nmax表示最大列号，height表示图像的高度，width表示图像的宽度。
进一步，所述步骤S312中搜索整条边界包括如下步骤：
S3121：新建与矩阵ImgIniSeg大小相同的零矩阵ImgBoundary。；
S3122：初始化最小行号mmin与最大行号mmax，以及最小列号nmin和最大列号nmax：
mmin＝mmax＝k＝m，nmin＝nmax＝p＝n，其中，m表示边界起点 ImgIniSeg(m,n)的行号，n表示边界起点ImgIniSeg(m,n)的列号，k表示以 ImgIniSeg(m,n)为起点的边界的下一个边界点的行号的变量，p为以为起点的边界的下一个边界点的列号的变量；
S3123：初始化提取边界方向，判断边界起点ImgIniSeg(m,n)的4-领域内，ImgIniSeg(m,n)的右方点值是否为1，若是，则设置方向dir＝0，同时p＝p+1, nmax＝p；若否，则设置方向dir＝1，同时k＝l+1,mmax＝k；同时，设置ImgBoundary (k,p)＝1，另记录起点的位置(k',p')和方向dir'，k'＝k,p'＝p,dir'＝dir；
S3124：将S3123中的方向dir的值代入dir＝(dir+3)％4，即计算dir+3的和除以4的余数，当dir＝(dir+3)％4时flag＝0，并进入步骤S3126；否则进入 S3125；
S3125：当dir＝0时，若ImgIniSeg(k,p+1)＝1，则flag＝1,p＝p+1，同时，若p>nmax，则nmax＝p，进入步骤S3126；否则继续下一行。
当dir＝1时，若ImgIniSeg(k+1,p)＝1，则flag＝1,k＝k+1，同时，若k>mmax，则mmax＝p，进入步骤S3126；否则继续下一行。
当dir＝2时，若ImgIniSeg(k,p-1)＝1，则flag＝1,p＝p-1，同时，若p<nmin，则nmin＝p，进入步骤S3126；否则继续下一行。
当dir＝3时，若ImgIniSeg(k-1,p)＝1，则flag＝1,k＝k-1，同时，若k<mmin，则mmin＝k，进入步骤S3126；
S3126：若flag＝1，则ImgBoundary(k,p)＝1，进入步骤S3127；若否，则 dir＝(dir+1)％4，进入步骤S3124；
S3127：判断条件(k＝k')&&(p＝p')&&(dir＝dir')是否满足，若满足，则表示该条边界搜索结束；若否，则进入步骤S311。
进一步，所述步骤S4中的提取的边界矩阵进行填充包括如下步骤：
S401：初始化矩阵ImgBoundaryFill，所述矩阵ImgBoundaryFill是与矩阵ImgBoundary大小相同的零矩阵；
S402：从上到下，从左到右的顺序扫描矩阵ImgBoundary，若该点的值为 0，且位于该点上方的点的值以及位于该点左方的点的值均不为0，则将该点视为土壤，同时将该点的值置为2；
S403：重复步骤S402，直至矩阵ImgBoundary遍历完成；
S404：从上到下，从左到右的顺序扫描矩阵ImgBoundary，寻找点 ImgBoundary(m，n)＝2且紧邻该点的右方或下方的点的值为0的点，将所述 ImgBoundary(m，n)置为0，同时，k＝m,p＝n,并将(k,p)压入栈stack；
S405：按右下左上的顺序在矩阵ImgBoundary中搜索一个在点(k,p)的4- 邻域内值为2的点，若搜索到新点，则停止下一个方向的搜索，更新(k,p)为新点的位置，ImgBoundary(k,p)＝0，将(k,p)压入栈stack；若未搜索到，从栈中出栈一个元素赋值给(k,p)；
S406：判断(k＝m)&&(p＝n)是否成立，若成立，则进入步骤S407；若否，则进入步骤S405；
S407：判断m＝height,n＝width是否成立，若成立，将矩阵ImgBoundaryFill 中与边界矩阵ImgBoundary中值不为0的点对应位置的点置为1，并结束，获得二值矩阵ImgBoundaryFill，即为二值图像III；若否，则进入步骤S404。
本发明的有益技术效果：本发明充分考虑到紫色土彩色图像的土壤区域在 Lab颜色空间L分量、a分量、b分量独立性高更适合做目标显著性特征表达，将紫色土彩色图像转换为Lab颜色空间，并在Lab颜色空间进行灰度变换，自适应聚类中心的密度峰值算法，克服了DPC算法中人工选择聚类中心容易出现失误的技术问题，将紫色土图像从含有背景区域的彩色图像中快速、准确、完整分割提取出来。
附图说明
下面结合附图和实施例对本发明作进一步描述：
图1为本发明的流程图。
图2为本发明中心决策度量CDMi示意图。
图3为本发明涉及传统的密度峰值聚类算法(DPC)决策图。
图4为本发明的紫色土原始图片。
图5为没有转换的原图的分割结果。
图6为CorrC2G方法转换后的图像分割结果。
图7为SeparableC2G方法转换后图像分割结果。
图8为不同算法对常规样本((No.5组))的分割结果。
图9不同算法对鲁棒样本(No.4组)的分割结果。
具体实施方式
以下结合说明书附图对本发明做出进一步的说明：
本发明提供的一种基于自适应密度峰值聚类的紫色土图像分割提取方法，
其特征在于：包括以下步骤：
S1：对含有紫色土区域的紫色土彩色图像进行可分离性灰度变换，获得灰度图像Ⅰ；
S2：利用自适应密度峰值聚类算法对灰度图像Ⅰ进行初步分割，获得初步分割后的二值图像Ⅱ；
S3：对二值图像Ⅱ进行边界提取处理，获得紫色土壤区域的边界矩阵；
S4：对提取的边界矩阵进行填充，获得二值图像Ⅲ；
S5：求出二值图像Ⅲ与所述含有紫色土区域的彩色图像的哈达玛积，得到只包含有紫色土图像的分割图像。
土壤“心土”未受到日晒雨淋等自然因素的影响，能充分显现土壤的颜色、结构等自然属性，是人工识别土壤信息的特征。由于本专利所做工作是机器视觉识别土壤种类的前期基础工作，因此，使用紫色土“心土”图像为试验样本。
试验图像采集方法：在野外，用铁锹锹出耕层10～20cm的紫色土，拍摄时需保证紫色土“心土”自然断口覆盖摄像头中心位置，并占全图至少50％的比例，获得紫色土“心土”自然断口图像。本专利算法在进行试验时，为了减小计算规模，将原始图像压缩为300×300像素大小。
首先对原始图像进行可分离性灰度变换获得灰度图像，然后利用自适应密度峰值聚类初步分割提取紫色土土壤目标，最后进行分割后处理得到精确的紫色土土壤区域图像。
本发明充分考虑到紫色土彩色图像的土壤区域在Lab颜色空间L分量、a分量、b分量独立性高更适合做目标显著性特征表达，将紫色土彩色图像转换为 Lab颜色空间，并在Lab颜色空间进行灰度变换，自适应聚类中心的密度峰值算法，克服了DPC算法中人工选择聚类中心容易出现失误的技术问题，将紫色土图像从含有背景区域的彩色图像中快速、准确、完整分割提取出来。
在本实施例中，所述步骤S1包括：
S101：将表征含有紫色土区域的紫色土彩色图像的矩阵ImgIn从RGB颜色空间转换到Lab颜色空间，获得矩阵ImgIn'，基于Lab颜色空间计算图像相似度矩阵Sx,y，以及矩阵ImgIn'中点ImgIn＇x高相似度点集合和低相似度点集合；
S102：在Lab颜色空间中建立灰度变换优化模型，采用离散搜索求解该优化模型，获得灰度矩阵ImgGray，即灰度图像Ⅰ。
由于Lab颜色空间的3个分量独立性高，更适合做目标显著性的特征表达。因此，在Lab颜色空间计算图像相似度矩阵，通过建立灰度变换优化模型获得能够改善聚类特征的可分离性灰度图像矩阵。
所述步骤S101包括：
S1011：在Lab颜色空间中，获取L、a、b分量的离散二维信息熵，所述L、 a、b分量的离散二维信息熵联立式子(1)至(3)确定，
Figure RE-GDA0002360084660000121
Figure RE-GDA0002360084660000122
c＝{L,a,b} (3)
其中，Hc表示单分量图像的离散二维熵，
Figure RE-GDA0002360084660000125
表示中间变量，
Figure RE-GDA0002360084660000123
表示单分量图像像素灰度值d与像素8邻域灰度均值e构成的二元组(d,e)出现的频数， eight表示图像的高度，width表示图像的宽度；
S1012：计算矩阵ImgIn'的相似度矩阵Sx,y，所述相似度矩阵Sx,y采用如下方法确定：
Figure RE-GDA0002360084660000124
其中，Sx,y表示矩阵ImgIn'的相似度矩阵，x＝1,...,width*height， y＝1,...,width*height，Hc表示单分量图像的离散二维熵，ImgIn＇c,x表示矩阵ImgIn'中点ImgIn＇x的c分量的值，ImgIn＇c,y表示矩阵ImgIn'中点ImgIn＇y的c分量的值；
Sx,y∈[0,1]，当Sx,y趋近于1时，表示点Ix和Iy相似度越大，越接近；反之，当 Sx,y趋近于0时，表示点Ix和Iy相似度越小，越远离；
S1013：计算矩阵ImgIn'中点ImgIn＇x高相似度点集合和低相似度点集合，
所述高相似度点集合采用如下方法确定：
set_highx＝{ImgIn＇y|Sx,y≥meanSx} (5)
其中，set_highx表示点ImgIn＇x的高相似度点集合，ImgIn＇y表示矩阵ImgIn'中除ImgIn＇x外的其他点，Sx,y表示矩阵ImgIn'的相似度矩阵，meanSx表示点ImgIn＇x与所有点间相似度的均值，
所述低相似度点集合采用如下方法确定：
set_lowx＝{ImgIn＇y|Sx,y＜meanSx} (6)
其中，set_lowx表示点ImgIn＇x的低相似度点集合，ImgIn＇y表示矩阵ImgIn'中除ImgIn＇x外的其他点，Sx,y表示矩阵ImgIn'的相似度矩阵，meanSx表示点ImgIn＇x与所有点间相似度的均值。
所述步骤S102包括：
S1021：建立优化模型计算Lab颜色空间亮度分量L的权重wL，颜色分量 a的权重wa,颜色分量b的权重wb,所述优化模型如下：
Figure RE-GDA0002360084660000131
其中，wL表示亮度分量L的权重，wa表示颜色分量a的权重，wb表示颜色分量b的权重，mheight表示图像的高度，nwidth表示图像的宽度,λ表示调节参数λ∈(0，G_hx和G_lx表示中间变量，
Figure RE-GDA0002360084660000132
Figure RE-GDA0002360084660000133
其中，ImgIn＇y表示矩阵ImgIn'中除ImgIn＇x外的其他点，set_highx表示点ImgIn＇x的高相似度点集合，set_lowx表示点ImgIn＇x的低相似度点集合，Sx,y表示矩阵ImgIn'的相似度矩阵， ImgGrayx表示任意点ImgIn＇x灰度变换后的值，ImgGrayy表示任意点ImgIn＇y灰度变换后的值，L表示Lab颜色空间中的亮度分量，a表示Lab颜色空间中的颜色分量，b表示Lab颜色空间中的亮度分量，wc表示单分量的权重，ImgIn'c,x表示点ImgIn＇x的某单分量值；
S1022：离散搜索求解优化模型(7)，获得紫色土灰度矩阵ImgGray，对所述紫色土灰度矩阵ImgGray采用如下方法获得：
a：求解临时矩阵ImgGrayTemp，
ImgGrayTempx＝wL*ImgIn＇L,x+wa*ImgIn＇a,x+wb*ImgIn＇b,x (8)
其中，ImgGrayTempx表示临时矩阵ImgGrayTemp中点x的值，wL表示亮度分量L的权重，wa表示颜色分量a的权重，wb表示颜色分量b的权重，ImgIn'L,x表示点ImgIn＇x的L分量值，ImgIn'a,x表示点ImgIn＇x的a分量值，ImgIn'b,x表示点ImgIn＇x的b分量值；
b：对矩阵ImgGrayTemp中元素向下取整，获得紫色土灰度矩阵ImgGray,即为灰度图像I。
对彩色图进行灰度变换是为了便于计算密度峰值，从而增大土壤和背景间的可分离性，便于后续图像分割。
在本实施例中，所述步骤S2包括：
S201：计算点ImgGraySi的密度ρi和点ImgGraySi的距离δi，其中，ImgGraySi表示数据集ImgGrayS中第i个点；
S202：自助确定聚类中心；
S203：利用自适应密度峰值聚类算法分割紫色土灰度矩阵ImgGray，获得二值紫色土初分割矩阵ImgIniSeg，即二值图像Ⅱ。
原始密度峰值聚类(density peaks clustering,DPC)是一种简洁、高效的聚类算法，被广泛应用于图像分割。DPC算法基于2个假设：聚类中心是密度ρ相对较大的数据点；聚类中心与比其密度大的其他数据点的距离δ相对较大。人工辅助选择ρ和δ都相对较大的点作为聚类中心，然后按照密度降序的顺序将每个剩余点分配给比其密度大的最近邻所属的类，得到任意形状的类。但其也有缺陷：根据经验确定参数截断距离dc，若参数不合适，则聚类效果较差；人工辅助选择聚类中心容易失误且不能满足野外自适应分割的需求。因此，本发明改进密度峰值聚类，提出自适应聚类中心的密度峰值聚类算法，实现紫色土的自适应分割。
所述步骤S201包括：
S2011：统计紫色土灰度矩阵ImgGray的一维直方图，直方图能够描述一幅图像中颜色的全局分布，而且容易理解和实现，所以在本实施例中都是基于紫色土灰度矩阵ImgGray的一维直方图，图像相似度计算都是使用它的，将ImgGray的灰度级作为待聚类数据集ImgGrayS＝{ImgGrayS1,…,ImgGraySi…,ImgGraySN},包含N 个数据，N＝max(ImgGray)-min(ImgGray),min(ImgGray)表示频数不为0的最小灰度级， max(ImgGray)表示频数不为0的最大灰度级；
S2012：计算点ImgGraySi的密度ρi，并记录ρi序排序时对应的原序号orderρi，其中密度ρi采用如下方法确定，
Figure RE-GDA0002360084660000151
其中，ρi表示点ImgGraySi的密度，N表示聚类数据集ImgGrayS的元素个数，frequencyj表示ImgGraySj的频数，dij为ImgGraySi与ImgGraySj间的欧式距离；
S2013：计算点ImgGraySi的距离δi，记录比ImgGraySi密度大的最近相邻点 ni，并将ImgGraySj赋值给ni，即ni＝ImgGraySj，对于密度最大的点，其ni＝0，其中，距离δi采用如下方法确定：
Figure RE-GDA0002360084660000152
其中，δi表示ImgGraySi的距离，dij为ImgGraySi与ImgGraySj间的欧式距离，δj表示非点ImgGraySi的点ImgGraySj的距离。
密度ρi是全局密度，距离越近的点对其贡献值越大，消除了DPC算法中数截断距离dc的影响。
所述步骤S202包括：
S2021：对ρi和δi分别进行min-max归一化处理，获得ρi'和δi'；定义任意点ImgGraySi到以ρi'为横轴，以δi'为纵轴构成的二维矩阵的左下角原点的线段作为三角形的斜边Hi，定义ImgGraySi到矩阵副对角线的垂线段作为三角形的对边 OSi，定义垂足点到原点的线段作为三角形的邻边ASi，
OSi＝sinβi×ISi (12)
ASi＝cosβi×ISi (13)
其中，δi'表示进行了min-max归一化处理ImgGrayiS的距离，ρi'表示进行了 min-max归一化处理点ImgGraySi的密度，Hi表示任意点ImgGraySi到以ρi'为横轴，以δi'为纵轴构成的二维矩阵的左下角原点的线段作为三角形的斜边，OSi表示ImgGraySi到矩阵副对角线的垂线段作为三角形的对边，ASi表示垂足点到原点的线段作为三角形的邻边，βi表示Hi与副对角线的夹角βi＝|αi-45|，Hi与横轴的夹角αi＝arctan(δi/ρi)/π*180；
S2022：计算中心决策度量CDMi，其中，中心决策度量CDMi采用如下方法确定，
CDMi＝Hi-OSi (14)
其中，Hi表示三角形的斜边，OSi表示三角形的对边；
S2023：对CDMi采用最大类间方差法确定分类阈值T，当CDMi≥T时，对应灰度级点即聚类中心点；在本实施例中，若聚类中心点数不足2个，将CDMi降序排序，取前面2个CDMi对应的灰度级点即为聚类中心点。
S2024：初始化每个点的类标记Ci，若ImgGraySi是聚类中心点，ci＝k， k表示第k个类，否则ci＝-1；
S2025：按orderρi的顺序将非聚类中心点归类到比其高密度点的最近邻点ni所属的类中：若ci＝-1，则ci＝cni，Ci表示数据集ImgGrayS中第i个点 ImgGraySi的类标记，ni表示比ImgGrayiS密度大的最近相邻点，其计算见步骤 S2013，因此，cni表示比ImgGraySi密度大的最近相邻点的类标记；
S2026：创建一个与矩阵ImgGray大小相同的零矩阵ImgIniSeg，遍历ImgGray全部元素，若ImgGraym,n＝ImgGraysi，则ImgIniSegm,n＝ci，m为行号，n为列号；
S2027：在矩阵ImgGray中心点附近随机抽取7个3×3大小子阵，剔除均值最大和最小的子阵，计算剩余子阵全部元素的均值μ作为土壤中心土块均值；计算每一类的聚类中心与土壤中心土块均值μ的距离，距离最小的那一类，将矩阵ImgIniSeg对应位置元素置1，否则置0，获得初分隔矩阵ImgIniSeg，即为初分割二值图像II。
在传统密度峰值聚类算法中，聚类中心点位于图3ρ-δ决策图右上角虚框内。根据CDMi定义，CDMi较大等效于ρ-δ决策图右上角虚框内点，则可将ρ-δ决策图右上角的定性决策转化为CDM的定量决策确定聚类中心点。即对CDM进行最大类间方差(OTSU)分类聚类，获得分类阈值T，CDMi≥T对应灰度级点即为聚类中心点。若聚类中心点数不足2个，将CDMi降序排序，取前面2个CDMi对应的灰度级点即为聚类中心点。因此，获得自适应聚类中心的密度峰值聚类算法(density peaks clustering with adaptive clustering centers,ACCDPC)。图2所示为中心决策量CDMi示意图。
在本实施例中，所述步骤S3中的边界提取包括以下步骤：
S311：对初分割矩阵ImgIniSeg中，搜索边界起点：扫描初分隔矩阵 ImgIniSeg，寻找边界起点ImgIniSeg(m,n)，所述边界起点ImgIniSeg(m,n)需同时满足：ImgIniSeg(m,n)＝1且在ImgIniSeg(m,n)的4-领域内，ImgIniSeg(m,n)的右方或下方的点至少有一个的值为1；
S312：搜索整条边界并记录该条边界最小行号mmin与最大行号mmax，以及该条边界最小列号nmin和最大列号nmax；
S313：若最小行号mmin与最大行号mmax，以及最小列号nmin和最大列号nmax满足以下式子，若满足，则结束；若不满足，将矩阵ImgIniSeg和ImgBoundary 中行号在mmin和mmax之间且列号在号nmin和nmax之间值为1的点置为0，并进入步骤S311；
Figure RE-GDA0002360084660000181
其中，mmin表示最小行号，mmax表示最大行号，nmin表示最小列号，nmax表示最大列号，height表示图像的高度，width表示图像的宽度。
所述步骤S312中搜索整条边界包括如下步骤：
S3121：新建与矩阵ImgIniSeg大小相同的零矩阵ImgBoundary。；
S3122：初始化最小行号mmin与最大行号mmax，以及最小列号nmin和最大列号nmax：
mmin＝mmax＝k＝m，nmin＝nmax＝p＝n，其中，m表示边界起点 ImgIniSeg(m,n)的行号，n表示边界起点ImgIniSeg(m,n)的列号，k表示以ImgIniSeg(m,n)为起点的边界的下一个边界点的行号的变量，p为以为起点的边界的下一个边界点的列号的变量；
S3123：初始化提取边界方向，判断边界起点ImgIniSeg(m,n)的4-领域内，ImgIniSeg(m,n)的右方点值是否为1，若是，则设置方向dir＝0，同时p＝p+1, nmax＝p；若否，则设置方向dir＝1，同时k＝l+1,mmax＝k；同时，设置ImgBoundary (k,p)＝1，另记录起点的位置(k',p')和方向dir'，k'＝k,p'＝p,dir'＝dir；
S3124：将S3123中的方向dir的值代入dir＝(dir+3)％4，即计算dir+3的和除以4的余数，当dir＝(dir+3)％4时flag＝0，并进入步骤S3126；否则进入 S3125；
S3125：当dir＝0时，若ImgIniSeg(k,p+1)＝1，则flag＝1,p＝p+1，同时，若p>nmax，则nmax＝p，进入步骤S3126；否则继续下一行。
当dir＝1时，若ImgIniSeg(k+1,p)＝1，则flag＝1,k＝k+1，同时，若k>mmax，则mmax＝p，进入步骤S3126；否则继续下一行。
当dir＝2时，若ImgIniSeg(k,p-1)＝1，则flag＝1,p＝p-1，同时，若p<nmin，则nmin＝p，进入步骤S3126；否则继续下一行。
当dir＝3时，若ImgIniSeg(k-1,p)＝1，则flag＝1,k＝k-1，同时，若k<mmin，则mmin＝k，进入步骤S3126；
S3126：若flag＝1，则ImgBoundary(k,p)＝1，进入步骤S3127；若否，则 dir＝(dir+1)％4，进入步骤S3124；
S3127：判断条件(k＝k')&&(p＝p')&&(dir＝dir')是否满足，若满足，则表示该条边界搜索结束；若否，则进入步骤S311。
在本实施例中，所述步骤S4中的提取的边界矩阵进行填充包括如下步骤：
S401：初始化矩阵ImgBoundaryFill，所述矩阵ImgBoundaryFill是与矩阵ImgBoundary大小相同的零矩阵；
S402：从上到下，从左到右的顺序扫描矩阵ImgBoundary，若该点的值为 0，且位于该点上方的点的值以及位于该点左方的点的值均不为0，则将该点视为土壤，同时将该点的值置为2；
S403：重复步骤S402，直至矩阵ImgBoundary遍历完成；
S404：从上到下，从左到右的顺序扫描矩阵ImgBoundary，寻找点 ImgBoundary(m，n)＝2且紧邻该点的右方或下方的点的值为0的点，将所述 ImgBoundary(m，n)置为0，同时，k＝m,p＝n,并将(k,p)压入栈stack；
S405：按右下左上的顺序在矩阵ImgBoundary中搜索一个在点(k,p)的4- 邻域内值为2的点，若搜索到新点，则停止下一个方向的搜索，更新(k,p)为新点的位置，ImgBoundary(k,p)＝0，将(k,p)压入栈stack；若未搜索到，从栈中出栈一个元素赋值给(k,p)；
S406：判断(k＝m)&&(p＝n)是否成立，若成立，则进入步骤S407；若否，则进入步骤S405；
S407：判断m＝height,n＝width是否成立，若成立，将矩阵ImgBoundaryFill 中与边界矩阵ImgBoundary中值不为0的点对应位置的点置为1，并结束，获得二值矩阵ImgBoundaryFill，即为二值图像III；若否，则进入步骤S404。
由于野外拍摄场景不可控，初分割图像土壤区域因光照、阴影、杂质等影响而包含内部空洞，背景区域可能包含同种土壤的离散小土块。因此，提出土壤边界提取及区域填充的后处理算法以获得完整的土壤区域图像，即本实施例中，所述步骤S3和S4。
为了验证本专利提供的基于自适应密度峰值聚类的紫色土图像分割提取方法的有效性，设计如下对比试验，并利用F1测度作为图像分割精度评价指标，其中F1＝2×P×R/(P+R)，其中P＝TP/(TP+FP)，R＝TP/(TP+FN)，TP为真阳性像素个数，FP为假阳性像素个数，FN为假阴性像素个数，TN为真阴性像素个数；以 (PS)手工分割的紫色土区域图像为参照标准。F1测度越大，图像分割精度越大。
对比试验如下：
(1)灰度变换对比试验：为验证可分离性灰度变换SeparableC2G算法的有效性，分别对其原始RGB彩色图像、彩色图像CorrC2G灰度变换后灰度图像、彩色图像经本专利SeparableC2G灰度变换后的灰度图像进行K-means算法 (k＝3)分割做对比试验。
(2)分割效果对比试验：为验证不同的分割算法对RGB彩色图像经 SeparableC2G灰度变换后的灰度图像进行分割的有效性，采用FRFCM算法、 H阈值分割算法、DPC算法、基于数据场的DFDPC算法共4种算法与本专利 ACCDPC聚类算法进行分割对比试验获得紫色土图像初分割结果；然后再用本专利后处理算法提取ACCDPC算法初分割结果的土壤区域边界并填充边界内土壤区域获得紫色土图像精确分割结果。
以上试验均在MatLab R2014a软件平台上和如下硬件配置上进行：Intel(R) Xeon(R)CPU E5-2687W v2@3.40GHz(2CPU)，内存64GB，显卡NVIDIA Quadro K5000的图形工作站，Windows 7专业版。
其中，试验对象为，从采集到的紫色土彩色图像中，随机抽取光照正常、周围无撒落“心土”粘连和无阴影遮蔽图像60张，构成20组常规样本；随机抽取背景中表土与“心土”颜色接近、周围有撒落“心土”粘连和阴影遮蔽图像60张，构成20组鲁棒样本。
随机选取(No.5组)常规样本和(No.4组)鲁棒样本显示本专利算法分割的图像结果。
(1)灰度变换效果对比试验结果与分析：
其对常规样本(No.5组)的图像进行灰度变换对比试验，分割结果如图4- 图7所示。
图4-图7所示：原始彩色图像或者CorrC2G图像的K-means算法分割结果均为离散的土块或颗粒且包含大量背景块，尤其是不能排除与目标土壤区域颜色非常接近的表土的影响，无法获得较完整紫色土土壤区域(图5、6)，原因在于原始彩色图像或者一般的灰度变换图像的土壤与背景间的差异性未达到能够将土壤区域与背景区域分割开来的程度。而SeparableC2G灰度图像来自于求解基于类间方差最大化类内方差最小化准则的优化模型，它提高了土壤区域与背景区域间的可分离性。先进行本专利SeparableC2G灰度变换后再K-means 分割可以获得紫色土土壤区域，但分割结果中包含许多空洞，对包含阴影的土壤图像分割不够完整(图7)。
(2)分割效果对比试验结果与分析：
对常规样本(No.5组)和鲁棒样本(No.4组)的图像分割结果分别如图8 和图9所示，全部20组常规样本和20组鲁棒样本的统计结果分别见表1和表 2。
如图8显示：对于常规样本，本专利ACC-DPC算法与4种对比算法基本都能分割出土壤区域。本专利ACC-DPC算法与DPC算法、DFDPC算法的分割结果较好且较为接近(图8d、图8e、图8f)，而FRFCM算法、H阈值算法受土壤表面凹凸不平造成不可避免的小而多的阴影块、表土和杂质等干扰因素影响较大，分割结果中土壤区域的空洞与背景区域的离散土块更多(图8b、图 8c)。ACCDPC算法加后处理算法，进一步去除背景区域的离散小土块和土壤区域的内部空洞(图8g)，与手工分割结果(图8h)接近。
如图9显示：对于鲁棒性样本，4种对比算法都出现不同程度的土壤分割不完整或背景区域粘连大块“心土”或“表土”，本专利ACC-DPC算法分割效果最好。FRFCM算法、DFDPC算法受阴影影响较大，导致土壤分割结果不完整，尤其是第三组图像这种土壤与背景都包含大片连续阴影图像表现更为明显(图9b、图9e)。H阈值算法对于鲁棒样本的适应能力最差，在对第一组图像的分割结果中该算法存在的土壤内部空洞问题最多，在对第三组图像的分割结果中该算法存在的将背景区域粘连大块“心土”或“表土”错分为土壤区域问题最严重(图9c)。DPC算法也对三组图像背景区域粘连大块“心土”或“表土”错分较多(图9d)。本专利ACCDPC算法的分割效果较好，背景错分为土壤的部分更少(图9f)。另外，DPC算法需要人工选择聚类中心，易出现失误，而ACC-DPC算法是自适应的，分割性能较为稳定。
表1展示了20组常规样本的每一组的平均分割精度和平均运行时间。从最后2行各算法对20组共60张常规样本的数据可以看出：(1)对于平均分割精度，FRFCM算法、H阈值算法、DPC算法、DFDPC算法、本专利ACCDPC 算法的值分别为81.91％、85.49％、90.29％、92.60％和93.45％，ACCDPC算法的平均分割精度最高，并且均方差最小，为2.17％，算法分割性能较稳定。另外，ACCDPC分割算法加上本专利边界提取和区域填充的后处理算法的平均分割精度为96.30％，进一步提升了图像分割精度。(2)对于平均运行时间，FRFCM 算法、H阈值算法、DPC算法、DFDPC算法、本专利ACCDPC算法的值分别为1.4s、0.42s、2.45s(含人工选择聚类中心耗时2.21s)、0.38s、0.26s，ACCDPC 算法耗时最少，而DPC算法耗时最久，原因在于该算法的人工选择聚类中心过程极大地增加了算法耗时，且不稳定。ACCDPC分割算法加上本专利后处理算法的平均运行时间为0.36s，可以看出后处理算法速度较快，耗时约0.1s。
表2展示了20组鲁棒样本的每一组的平均分割精度和平均运行时间。从最后2行各算法对20组共60张鲁棒样本的数据可以看出：与表2一样，本专利 ACCDPC算法在平均分割精度和运行时间两方面均优于4种对比算法，且与4 种对比算法在平均分割精度方面进一步拉开了差距。ACCDPC算法比FRFCM 算法、H阈值算法、DPC算法和DFDPC算法分别高出11.75％、5.2％、12.47％和3.09％，表明ACCDPC算法在分割背景更为复杂、光照影响更大的鲁棒性实验图像时优势更为明显。后处理算法进一步提升了本专利算法的分割精度。
表1 20组常规样本的统计值
Table 1 Statistical value of 20 groups of normal samples
Figure RE-GDA0002360084660000231
Figure RE-GDA0002360084660000241
注：*后的数字表示在DPC算法中，手工辅助选择聚类中心部分耗时。
表2 20组鲁棒样本的统计值
Table 2 Statistical value of 20 groups of robustness samples
注：*后的数字表示在DPC算法中，手工辅助选择聚类中心部分耗时。
表2最后两行数据与表1最后两行数据相比，FRFCM算法、H阈值算法、 DPC算法、DFDPC算法、本专利ACCDPC算法在平均分割精度方面都有所减小，说明非均匀光照带来的阴影、颜色相近的“表土”或“心土”等情况对图像分割结果都有一定影响，使得土壤区域出现被阴影覆盖的土壤部分缺失或者多出部分“表土”，因此分割精度下降。但从整体情况来看，ACCDPC算法的平均分割精度依旧是最高的，对复杂背景情况适应能力最强。
本专利提出了一种野外紫色土彩色图像分割方法，该方法首先对原始图像进行可分离性灰度变换获得灰度图像，然后利用本专利自适应密度峰值聚类初步分割提取土壤目标，最后利用土壤区域边界提取与区域填充的后处理算法获得精确的紫色土土壤区域图像。20组常规样本和20组鲁棒样本结果表明：
1)本专利算法的平均分割精度最高。本专利ACCDPC算法的初分割平均分割精度为93.45％(常规样本)和87.40％(鲁棒样本)，比FRFCM算法、H 阈值算法、DPC算法、DFDPC算法分别高出11.54％(常规样本)和11.75％(鲁棒样本)、7.96％(常规样本)和5.2％(鲁棒样本)、3.16％(常规样本)和 12.47％(鲁棒样本)、0.85％(常规样本)和3.09％(鲁棒样本)。由鲁棒样本结果可知，阴影、颜色相近的“表土”或“心土”等情况对图像分割结果有一定影响，但ACCDPC算法的平均分割精度依旧是最高的，说明本专利ACCDPC 算法对复杂背景情况适应能力最强。本专利后处理算法能够有效去除背景区域的离散小土块和土壤区域的内部空洞，进一步提高了平均分割精度，分别为 96.30％(常规样本)和91.63％(鲁棒样本)。
2)本专利算法的平均耗时最短。常规样本结果中，FRFCM算法、H阈值算法、DPC算法、DFDPC算法、本专利ACCDPC算法初分割平均耗时分别为 1.4s、0.42s、2.45s(含人工选择聚类中心耗时2.21s)、0.38s、0.26s，ACCDPC 算法耗时最少。后处理算法的时间花销约0.1s，后处理算法是快速有效的。鲁棒样本的耗时情况与常规样本基本一致。
最后说明的是，以上实施例仅用以说明本发明的技术方案而非限制，尽管参照较佳实施例对本发明进行了详细说明，本领域的普通技术人员应当理解，可以对本发明的技术方案进行修改或者等同替换，而不脱离本发明技术方案的宗旨和范围，其均应涵盖在本发明的权利要求范围当中。

Claims (10)

1.一种基于自适应密度峰值聚类的紫色土图像分割提取方法，其特征在于：
包括以下步骤：
S1：对含有紫色土区域的紫色土彩色图像进行可分离性灰度变换，获得灰度图像Ⅰ；
S2：利用自适应密度峰值聚类算法对灰度图像Ⅰ进行初步分割，获得初步分割后的二值图像Ⅱ；
S3：对二值图像Ⅱ进行边界提取处理，获得紫色土壤区域的边界矩阵；
S4：对提取的边界矩阵进行填充，获得二值图像Ⅲ；
S5：求出二值图像Ⅲ与所述含有紫色土区域的彩色图像的哈达玛积，得到只包含有紫色土图像的分割图像。
2.根据权利要求1所述基于自适应密度峰值聚类的紫色土图像分割提取方法，其特征在于：所述步骤S1包括：
S101：将表征含有紫色土区域的紫色土彩色图像的矩阵ImgIn从RGB颜色空间转换到Lab颜色空间，获得矩阵ImgIn'，基于Lab颜色空间计算图像相似度矩阵Sx,y，以及矩阵ImgIn'中点ImgIn＇x高相似度点集合和低相似度点集合；
S102：在Lab颜色空间中建立灰度变换优化模型，采用离散搜索求解该优化模型，获得灰度矩阵ImgGray，即灰度图像Ⅰ。
3.根据权利要求2所述基于自适应密度峰值聚类的紫色土图像分割提取方法，其特征在于：所述步骤S101包括：
S1011：在Lab颜色空间中，获取L、a、b分量的离散二维信息熵，所述L、a、b分量的离散二维信息熵联立式子(1)至(3)确定，
Figure FDA0002238322710000011
Figure FDA0002238322710000021
c＝{L,a,b} (3)
其中，Hc表示单分量图像的离散二维熵，p(c) de表示中间变量，
Figure FDA0002238322710000022
表示单分量图像像素灰度值d与像素8邻域灰度均值e构成的二元组(d,e)出现的频数，eight表示图像的高度，width表示图像的宽度；
S1012：计算矩阵ImgIn'的相似度矩阵Sx,y，所述相似度矩阵Sx,y采用如下方法确定：
Figure FDA0002238322710000023
其中，Sx,y表示矩阵Img'I的相似度矩阵，x＝1,...,width*height，y＝1,...,width*height，Hc表示单分量图像的离散二维熵，ImgIn＇c,x表示矩阵ImgIn'中点ImgIn＇x的c分量的值，ImgIn＇c,y表示矩阵ImgIn'中点ImgIn＇y的c分量的值；
S1013：计算矩阵ImgIn'中点ImgIn＇x高相似度点集合和低相似度点集合，
所述高相似度点集合采用如下方法确定：
set-highx＝{ImgIn＇y|Sx,y≥meanSx} (5)
其中，set-highx表示点ImgIn＇x的高相似度点集合，ImgIn＇y表示矩阵ImgIn'中除ImgIn＇x外的其他点，Sx,y表示矩阵ImgIn'的相似度矩阵，meanSx表示点ImgIn＇x与所有点间相似度的均值，
所述低相似度点集合采用如下方法确定：
set-lowx＝{ImgIn＇y|Sx,y＜meanSx} (6)
其中，set-lowx表示点ImgIn＇x的低相似度点集合，ImgIn＇y表示矩阵ImgIn'中除ImgIn＇x外的其他点，Sx,y表示矩阵ImgIn'的相似度矩阵，meanSx表示点ImgIn＇x与所有点间相似度的均值。
4.根据权利要求2所述基于自适应密度峰值聚类的紫色土图像分割提取方法，其特征在于：所述步骤S102包括：
S1021：建立灰度变换优化模型计算Lab颜色空间亮度分量L的权重wL，颜色分量a的权重wa,颜色分量b的权重wb,所述优化模型如下：
Figure FDA0002238322710000031
其中，wL表示亮度分量L的权重，wa表示颜色分量a的权重，wb表示颜色分量b的权重，height表示图像的高度，width表示图像的宽度,λ表示调节参数λ∈(，G_hx和G_lx表示中间变量，
Figure FDA0002238322710000033
其中，ImgIn＇y表示矩阵ImgIn'中除ImgIn＇x外的其他点，set-highx表示点ImgIn＇x的高相似度点集合，set-lowx表示点ImgIn＇x的低相似度点集合，Sx,y表示矩阵ImgIn'的相似度矩阵，ImgGrayx表示任意点ImgIn＇x灰度变换后的值，ImgGrayy表示任意点ImgIn＇y灰度变换后的值，L表示Lab颜色空间中的亮度分量，a表示Lab颜色空间中的颜色分量，b表示Lab颜色空间中的亮度分量，wc表示单分量的权重，ImgIn'c,x表示点ImgIn＇x的某单分量值；
S1022：离散搜索求解优化模型(7)，获得紫色土灰度矩阵ImgGray，对所述紫色土灰度矩阵ImgGray采用如下方法获得：
a：求解临时矩阵ImgGrayTemp，
ImgGrayTempx＝wL*ImgIn＇L,x+wa*ImgIn＇a,x+wb*ImgIn＇b,x (8)
其中，ImgGrayTempx表示临时矩阵ImgGrayTemp中点x的值，wL表示亮度分量L的权重，wa表示颜色分量a的权重，wb表示颜色分量b的权重，ImgIn'L,x表示点ImgIn＇x的L分量值，ImgIn'a,x表示点ImgIn＇x的a分量值，ImgIn'b,x表示点ImgIn＇x的b分量值；
b：对矩阵ImgGrayTemp中元素向下取整，获得紫色土灰度矩阵ImgGray,即为灰度图像I。
5.根据权利要求1所述基于自适应密度峰值聚类的紫色土图像分割提取方法，其特征在于：所述步骤S2包括：
S201：计算点ImgGraySi的密度ρi和点ImgGraySi的距离δi，其中，ImgGraySi表示数据集ImgGrayS中第i个点；
S202：自助确定聚类中心；
S203：利用自适应密度峰值聚类算法分割紫色土灰度矩阵ImgGray，获得二值紫色土初分割矩阵ImgIniSeg，即二值图像Ⅱ。
6.根据权利要求5所述基于自适应密度峰值聚类的紫色土图像分割提取方法，其特征在于：所述步骤S201包括：
S2011：统计紫色土灰度矩阵ImgGray的一维直方图，将ImgGray的灰度级作为待聚类数据集ImgGrayS＝{ImgGrayS1,…,ImgGraySi…,ImgGraySN},包含N个数据，N＝max(ImgGray)-min(ImgGray),min(ImgGray)表示频数不为0的最小灰度级，max(ImgGray)表示频数不为0的最大灰度级；
S2012：计算点ImgGraySi的密度ρi，并记录ρi序排序时对应的原序号orderρi，其中密度ρi采用如下方法确定，
Figure FDA0002238322710000041
其中，ρi表示点ImgGraySi的密度，N表示聚类数据集ImgGrayS的元素个数，frequencyj表示ImgGraySj的频数，dij为ImgGraySi与ImgGraySj间的欧式距离；
S2013：计算点ImgGraySi的距离δi，记录比ImgGraySi密度大的最近相邻点ni，并将ImgGraySj赋值给ni，即ni＝ImgGraySj，对于密度最大的点，其ni＝0，其中，距离δi采用如下方法确定：
其中，δi表示ImgGraySi的距离，dij为ImgGraySi与ImgGraySj间的欧式距离，δj表示非点ImgGraySi的点ImgGraySj的距离。
7.根据权利要求5所述基于自适应密度峰值聚类的紫色土图像分割提取方法，其特征在于：所述步骤S202包括：
S2021：对ρi和δi分别进行min-max归一化处理，获得ρi'和δi'；定义任意点ImgGraySi到以ρi'为横轴，以δi'为纵轴构成的二维矩阵的左下角原点的线段作为三角形的斜边Hi，定义ImgGraySi到矩阵副对角线的垂线段作为三角形的对边OSi，定义垂足点到原点的线段作为三角形的邻边ASi，
Figure FDA0002238322710000052
OSi＝sinβi×Hi (12)
ASi＝cosβi×Hi (13)
其中，δi'表示进行了min-max归一化处理后点ImgGraySi的距离，ρi'表示进行了min-max归一化处理后点ImgGraySi的密度，Hi表示任意点ImgGraySi到以ρi'为横轴，以δi'为纵轴构成的二维矩阵的左下角原点的线段作为三角形的斜边，OSi表示ImgGrayiS到矩阵副对角线的垂线段作为三角形的对边，ASi表示垂足点到原点的线段作为三角形的邻边，βi表示Hi与副对角线的夹角βi＝|αi-45|，Hi与横轴的夹角αi＝arctan(δi/ρi)/π*180；
S2022：计算中心决策度量CDMi，其中，中心决策度量CDMi采用如下方法确定，
CDMi＝Hi-OSi(14)
其中，Hi表示三角形的斜边，OSi表示三角形的对边；
S2023：对CDMi采用最大类间方差法确定分类阈值T，当CDMi≥T时，对应灰度级点即聚类中心点；
S2024：初始化每个点的类标记Ci，若ImgGraySi是聚类中心点，ci＝k，k表示第k个类，否则ci＝-1；
S2025：按orderρi的顺序将非聚类中心点归类到比其高密度点的最近邻点ni所属的类中：若ci＝-1，则ci＝cni，Ci表示数据集ImgGrayS中第i个点ImgGraySi的类标记，cni表示比ImgGraySi密度大的最近相邻点，其计算见步骤S3013；
S2026：创建一个与矩阵ImgGray大小相同的零矩阵ImgIniSeg，遍历ImgGray全部元素，若ImgGraym,n＝ImgGraysi，则ImgIniSegm,n＝ci，m为行号，n为列号；
S2027：在矩阵ImgGray中心点附近随机抽取7个3×3大小子阵，剔除均值最大和最小的子阵，计算剩余子阵全部元素的均值μ作为土壤中心土块均值；计算每一类的聚类中心与土壤中心土块均值μ的距离，距离最小的那一类，将矩阵ImgIniSeg对应位置元素置1，否则置0，获得初分隔矩阵ImgIniSeg，即为初分割二值图像II。
8.根据权利要求1所述基于自适应密度峰值聚类的紫色土图像分割提取方法，其特征在于：所述步骤S3中的边界提取包括以下步骤：
S311：对初分割矩阵ImgIniSe中，搜索边界起点：扫描初分隔矩阵ImgIniSeg，寻找边界起点ImgIniSeg(m,n)，所述边界起点ImgIniSeg(m,n)需同时满足：ImgIniSeg(m,n)＝1且在ImgIniSeg(m,n)的4-领域内，ImgIniSeg(m,n)的右方或下方的点至少有一个的值为1；
S312：搜索整条边界并记录该条边界最小行号mmin与最大行号mmax，以及该条边界最小列号nmin和最大列号nmax；
S313：判断最小行号mmin与最大行号mmax，以及最小列号nmin和最大列号nmax满足其中，mmin表示最小行号，mmax表示最大行号，nmin表示最小列号，nmax表示最大列号，height表示图像的高度，width表示图像的宽度，
若满足，则结束；若不满足，将矩阵ImgIniSeg和ImgBoundary中行号在mmin和mmax之间且列号在号nmin和nmax之间值为1的点置为0，并进入步骤S311。
9.根据权利要求8所述基于自适应密度峰值聚类的紫色土图像分割提取方法，其特征在于：所述步骤S312中搜索整条边界包括如下步骤：
S3121：新建与矩阵ImgIniSeg大小相同的零矩阵ImgBoundary；
S3122：初始化最小行号mmin与最大行号mmax，以及最小列号nmin和最大列号nmax：
mmin＝mmax＝k＝m，nmin＝nmax＝p＝n，其中，m表示边界起点ImgIniSeg(m,n)的行号，n表示边界起点ImgIniSeg(m,n)的列号，k表示以ImgIniSeg(m,n)为起点的边界的下一个边界点的行号的变量，p为以为起点的边界的下一个边界点的列号的变量；
S3123：初始化提取边界方向，判断边界起点ImgIniSeg(m,n)的4-领域内，ImgIniSeg(m,n)的右方点值是否为1，若是，则设置方向dir＝0，同时p＝p+1,nmax＝p；若否，则设置方向dir＝1，同时k＝l+1,mmax＝k；同时，设置ImgBoundary(k,p)＝1，另记录起点的位置(k',p')和方向dir'，k'＝k,p'＝p,dir'＝dir；
S3124：将S3123中的方向dir的值代入dir＝(dir+3)％4，即计算dir+3的和除以4的余数，当dir＝(dir+3)％4时flag＝0，并进入步骤S3126；否则进入S3125；
S3125：当dir＝0时，若ImgIniSeg(k,p+1)＝1，则flag＝1,p＝p+1，同时，若p>nmax，则nmax＝p，进入步骤S3126；否则继续下一行。
当dir＝1时，若ImgIniSeg(k+1,p)＝1，则flag＝1,k＝k+1，同时，若k>mmax，则mmax＝p，进入步骤S3126；否则继续下一行。
当dir＝2时，若ImgIniSeg(k,p-1)＝1，则flag＝1,p＝p-1，同时，若p<nmin，则nmin＝p，进入步骤S3126；否则继续下一行。
当dir＝3时，若ImgIniSeg(k-1,p)＝1，则flag＝1,k＝k-1，同时，若k<mmin，则mmin＝k，进入步骤S3126；
S3126：若flag＝1，则ImgBoundary(k,p)＝1，进入步骤S3127；若否，则dir＝(dir+1)％4，进入步骤S3124；
S3127：判断条件(k＝k')&&(p＝p')&&(dir＝dir')是否满足，若满足，则表示该条边界搜索结束；若否，则进入步骤S311。
10.根据权利要求1所述基于自适应密度峰值聚类的紫色土图像分割提取方法，其特征在于：所述步骤S4中的提取的边界矩阵进行填充包括如下步骤：
S401：初始化矩阵ImgBoundaryFill，所述矩阵ImgBoundaryFill是与矩阵ImgBoundary大小相同的零矩阵；
S402：从上到下，从左到右的顺序扫描矩阵ImgBoundary，若该点的值为0，且位于该点上方的点的值以及位于该点左方的点的值均不为0，则将该点视为土壤，同时将该点的值置为2；
S403：重复步骤S402，直至矩阵ImgBoundary遍历完成；
S404：从上到下，从左到右的顺序扫描矩阵ImgBoundary，寻找点ImgBoundary(m，n)＝2且紧邻该点的右方或下方的点的值为0的点，将所述ImgBoundary(m，n)置为0，同时，k＝m,p＝n,并将(k,p)压入栈stack；
S405：按右下左上的顺序在矩阵ImgBoundary中搜索一个在点(k,p)的4-邻域内值为2的点，若搜索到新点，则停止下一个方向的搜索，更新(k,p)为新点的位置，ImgBoundary(k,p)＝0，将(k,p)压入栈stack；若未搜索到，从栈中出栈一个元素赋值给(k,p)；
S406：判断(k＝m)&&(p＝n)是否成立，若成立，则进入步骤S407；若否，则进入步骤S405；
S407：判断m＝height,n＝width是否成立，若成立，将矩阵ImgBoundaryFill中与边界矩阵ImgBoundary中值不为0的点对应位置的点置为1，并结束，获得二值矩阵ImgBoundaryFill，即为二值图像III；若否，则进入步骤S404。
CN201910991121.XA 2019-10-18 2019-10-18 基于自适应密度峰值聚类的紫色土图像分割提取方法 Active CN110827306B (zh)

Priority Applications (1)

Application Number Priority Date Filing Date Title
CN201910991121.XA CN110827306B (zh) 2019-10-18 2019-10-18 基于自适应密度峰值聚类的紫色土图像分割提取方法

Applications Claiming Priority (1)

Application Number Priority Date Filing Date Title
CN201910991121.XA CN110827306B (zh) 2019-10-18 2019-10-18 基于自适应密度峰值聚类的紫色土图像分割提取方法

Publications (2)

Publication Number Publication Date
CN110827306A true CN110827306A (zh) 2020-02-21
CN110827306B CN110827306B (zh) 2022-06-21

Family

ID=69549468

Family Applications (1)

Application Number Title Priority Date Filing Date
CN201910991121.XA Active CN110827306B (zh) 2019-10-18 2019-10-18 基于自适应密度峰值聚类的紫色土图像分割提取方法

Country Status (1)

Country Link
CN (1) CN110827306B (zh)

Cited By (8)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN112348761A (zh) * 2020-11-30 2021-02-09 广州绿怡信息科技有限公司 设备外观图像亮度调节方法及装置
CN112434740A (zh) * 2020-11-26 2021-03-02 西北大学 一种基于深度学习的秦俑碎块分类方法
CN113205525A (zh) * 2021-03-23 2021-08-03 中国科学院南京土壤研究所 一种基于超绿指数及超像素的土壤图像分割提取方法
CN113808144A (zh) * 2021-09-22 2021-12-17 重庆师范大学 基于直方图拟合的密度峰值聚类土壤图像阴影检测方法
CN114994102A (zh) * 2022-08-04 2022-09-02 武汉钰品研生物科技有限公司 基于x光的食品异物无痕快速检测方法
CN116223515A (zh) * 2023-05-05 2023-06-06 成都中航华测科技有限公司 一种用于电路板测试过程的导电图形缺陷检测方法
CN116894851A (zh) * 2023-09-11 2023-10-17 吉林农业大学 一种土壤ct图像孔隙提取方法和系统
CN117237646A (zh) * 2023-11-15 2023-12-15 深圳市润海电子有限公司 基于图像分割的pet高温阻燃胶带瑕疵提取方法及系统

Citations (2)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN108961295A (zh) * 2018-07-27 2018-12-07 重庆师范大学 基于正态分布h阈值的紫色土图像分割提取方法
CN109035289A (zh) * 2018-07-27 2018-12-18 重庆师范大学 基于切比雪夫不等式h阈值的紫色土图像分割提取方法

Patent Citations (2)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN108961295A (zh) * 2018-07-27 2018-12-07 重庆师范大学 基于正态分布h阈值的紫色土图像分割提取方法
CN109035289A (zh) * 2018-07-27 2018-12-18 重庆师范大学 基于切比雪夫不等式h阈值的紫色土图像分割提取方法

Non-Patent Citations (4)

* Cited by examiner, † Cited by third party
Title
GONZALEZ R C等: "Digital Image Processing", 《ENGLEWOOD CLIFFS, NJ: PRENTICE HALL》 *
卢红阳等: "基于最大加权投影求解的彩色图像", 《自动化学报》 *
曾绍华等: "基于自适应密度峰值聚类的野外紫色土彩色图像分割", 《农业工程学》 *
贺付亮等: "基于视觉显著性和脉冲耦合神经网络", 《农业工程学报》 *

Cited By (12)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN112434740A (zh) * 2020-11-26 2021-03-02 西北大学 一种基于深度学习的秦俑碎块分类方法
CN112348761A (zh) * 2020-11-30 2021-02-09 广州绿怡信息科技有限公司 设备外观图像亮度调节方法及装置
CN112348761B (zh) * 2020-11-30 2024-04-26 广州绿怡信息科技有限公司 设备外观图像亮度调节方法及装置
CN113205525A (zh) * 2021-03-23 2021-08-03 中国科学院南京土壤研究所 一种基于超绿指数及超像素的土壤图像分割提取方法
CN113808144A (zh) * 2021-09-22 2021-12-17 重庆师范大学 基于直方图拟合的密度峰值聚类土壤图像阴影检测方法
CN113808144B (zh) * 2021-09-22 2024-03-29 重庆师范大学 基于直方图拟合的密度峰值聚类土壤图像阴影检测方法
CN114994102A (zh) * 2022-08-04 2022-09-02 武汉钰品研生物科技有限公司 基于x光的食品异物无痕快速检测方法
CN116223515A (zh) * 2023-05-05 2023-06-06 成都中航华测科技有限公司 一种用于电路板测试过程的导电图形缺陷检测方法
CN116894851A (zh) * 2023-09-11 2023-10-17 吉林农业大学 一种土壤ct图像孔隙提取方法和系统
CN116894851B (zh) * 2023-09-11 2023-11-14 吉林农业大学 一种土壤ct图像孔隙提取方法和系统
CN117237646A (zh) * 2023-11-15 2023-12-15 深圳市润海电子有限公司 基于图像分割的pet高温阻燃胶带瑕疵提取方法及系统
CN117237646B (zh) * 2023-11-15 2024-01-30 深圳市润海电子有限公司 基于图像分割的pet高温阻燃胶带瑕疵提取方法及系统

Also Published As

Publication number Publication date
CN110827306B (zh) 2022-06-21

Similar Documents

Publication Publication Date Title
CN110827306B (zh) 基于自适应密度峰值聚类的紫色土图像分割提取方法
CN108038476B (zh) 一种基于边缘检测与sift的人脸表情识别特征提取方法
Es-saady et al. Automatic recognition of plant leaves diseases based on serial combination of two SVM classifiers
CN108388874B (zh) 基于图像识别与级联分类器的对虾形态参数自动测量方法
CN110120042B (zh) 一种基于slic超像素和自动阈值分割的农作物图像病虫害区域提取方法
CN106340016B (zh) 一种基于细胞显微镜图像的dna定量分析方法
CN108629783B (zh) 基于图像特征密度峰值搜索的图像分割方法、系统及介质
CN108596195B (zh) 一种基于稀疏编码特征提取的场景识别方法
CN113536958B (zh) 导航路径提取方法、装置、农业机器人和存储介质
CN109426793A (zh) 一种图像行为识别方法、设备及计算机可读存储介质
CN110533648A (zh) 一种黑头识别处理方法及系统
CN110874835B (zh) 作物叶部病害抗性鉴定方法、系统、电子设备及存储介质
CN115690086A (zh) 一种基于对象的高分辨率遥感影像变化检测方法及系统
CN116612307A (zh) 一种基于迁移学习的茄科病害等级识别方法
CN115311689A (zh) 牛脸身份识别特征提取模型构建方法及牛脸身份识别方法
CN110084820B (zh) 基于改进fcm算法的紫色土图像自适应分割提取方法
CN111414958A (zh) 一种视觉词袋金字塔的多特征图像分类方法及系统
CN113657216A (zh) 基于形状特征的点云场景树木树冠及木质点分离方法
CN116596921B (zh) 一种焚烧炉渣分选方法及系统
US11847811B1 (en) Image segmentation method combined with superpixel and multi-scale hierarchical feature recognition
Rathore et al. CBISC: a novel approach for colon biopsy image segmentation and classification
CN109299295B (zh) 蓝印花布图像数据库搜索方法
CN109753912B (zh) 一种基于张量的多光谱掌纹匹配方法
CN115511838A (zh) 一种基于群智能优化的植物病害高精度识别方法
CN112686302B (zh) 图像特征点匹配方法、计算机设备及存储装置

Legal Events

Date Code Title Description
PB01 Publication
PB01 Publication
SE01 Entry into force of request for substantive examination
SE01 Entry into force of request for substantive examination
GR01 Patent grant
GR01 Patent grant