CN103065308B - 一种基于马尔科夫随机场的粘连物分割方法 - Google Patents

Abstract

本发明公开了质量监测技术领域，特别涉及一种基于马尔科夫随机场的粘连物分割方法。本发明将粘连物图像分割为正确分割区域和过分割区域；然后计算所述过分割区域和该所述过分割区域的八联通区域之间的紧密度，根据所述紧密度对所述过分割区域进行合并，得到粘连物的分割结果。本发明利用粘连物图像的面积和相对高度作为二维观测数据，结合图像的二维空间特征相关性，将过分割问题转化为求解图像区域标记场问题；引入分水岭变换解决过分割问题的基础上，依据邻域间的二维空间关系赋予不同势能，以抑制噪声影响，最后计算过分割区域与邻域的紧密度，选择紧密度最大的邻域并与之合并，解决过分割问题，从而实现粘连果蔬合格率的在线监测。

Description

一种基于马尔科夫随机场的粘连物分割方法

技术领域

本发明涉及质量监测技术领域，特别涉及一种基于马尔科夫随机场的粘连物分割方法。

背景技术

果蔬质量控制技术是指在不接触、不破坏鲜切果蔬的前提下，对果蔬的颜色、尺寸、形状等外部品质特征进行检测，其中的一个难点问题是果蔬之间的粘连导致的误判。

最近20年来，随着计算机技术的不断进步，机器视觉已经成为果蔬合格率在线监测的一种有效方法，但是大部分果蔬合格率在线监测研究建立在人为的将粘连果蔬分离的基础上，这是一个非常浪费时间的过程。因此，自动分割粘连图像是果蔬合格率在线监测的基础。

针对果蔬粘连问题，国内外学者提出了多种粘连图像分割方法，如：数学形态学方法，分水岭方法，傅里叶级数近似算法，基于曲率的算法，椭圆拟合方法等。

以上方法能在某些特定的场景中有好的应用，但是对于果蔬流水线上存在大量粘连的图像来说，有的由20个以上的果蔬粘连在一起，用上述方法很难得到满意的分割结果。因此，需要一种新的果蔬分级方法，以能够有效地对粘连果蔬进行分割，并进行合格率无损监测。

发明内容

(一)要解决的技术问题

本发明要解决的技术问题是：如何实现快速、准确、自动、无损地分割粘连果蔬。

(二)技术方案

为了解决上述技术问题，本发明提供一种基于马尔科夫随机场的粘连物分割方法，其特征是，该方法包括以下步骤：

S1：将粘连物图像分割为正确分割区域和过分割区域；

S2：计算所述过分割区域和该所述过分割区域的八联通区域之间的紧密度，根据所述紧密度对所述过分割区域进行合并，得到粘连物的分割结果。

所述步骤S1具体为：

S11：将粘连物图像转换为距离灰度图像；

S12：对所述距离灰度图像进行处理得到初始分水岭；

S13：通过所述初始分水岭得到所述距离灰度图像的分割区域；

S14：通过所述分割区域的二维观测特征将所述分割区域划分为正确分割区域或过分割区域。

所述步骤S12具体为：

通过分水岭法对所述距离灰度图像进行处理得到初始分水岭。

所述步骤S14具体为：

S141：获取所述分割区域的二维观测特征；

S142：根据所述二维观测特征将所述分割区域划分为标记区域；

S143：根据所述标记区域内的区域的观测特征数据计算所述标记区域内的区域的能量；

S144：重复执行步骤S141至S143，分别比较相邻两次得到的所述标记区域内的区域的能量之间的差；

S145：若所述能量之间的差小于设定阈值，则将该所述标记区域划分为正确分割区域或过分割区域。

所述标记区域内的区域的能量的计算公式为：

E＝E_L+E_F

其中：

E为标记区域内的区域的能量；

E_L为标记区域内的区域标记势能；

E_F为标记区域内的区域标记场能量。

所述的E_L计算公式为：

$E_L=\underset{S}{\mathrm{\Σ}}{\mathrm{\β}}_s\×\frac{k_r}{\underset{r\∈N_s}{\mathrm{\Σ}}{k}_r}\×\underset{r\∈N_s}{\mathrm{\Σ}}\mathrm{\δ}(L_s,L_r)$

其中：

S为所有标记区域的集合；

β_s为区域s的耦合系数；

k_r为区域s与区域r之间的紧密度；

N_s为区域s的所有邻域；

δ(L_s,L_r)为势团势能V_{s,r}(L_s,L_r)的符号函数。

所述势团势能的计算公式为：

$V_{\{s,r\}}(L_s,L_r)=\left\{\begin{array}{c}-{\mathrm{\β}}_s\×\frac{k_r}{\underset{r\∈N_s}{\mathrm{\Σ}}{k}_r},L_s=L_r\\ {\mathrm{\β}}_s\×\frac{k_r}{\underset{r\∈N_s}{\mathrm{\Σ}}{k}_r},L_s\≠L_r\end{array}\right.$

其中：

V_{s,r}(L_s,L_r)为标记为L_s和L_r的标记场之间的势团势能；

所述的E_F计算公式为：

$E_F=\underset{S,L_S=m}{\mathrm{\Σ}}\left\{\underset{k=1}{\overset{2}{\mathrm{\Σ}}}\right[\frac{{(f_s^k-{\mathrm{\μ}}_m^k)}^2}{2{\mathrm{\σ}}_m^{k2}}+\log \left(\sqrt{2\mathrm{\π}}{\mathrm{\σ}}_m^{\text{k}}\right)\left]\right\}$

其中：

S,L_S＝m为所有第m类区域的集合；

为第m类区域的第k维观测数据；

为标记区域内的区域的第k维观测数据的均值；

为标记区域内的区域的第k维观测数据的方差。

所述紧密度的计算公式为：

t＝αSim+βd+γp

其中：

t为紧密度；

α为灰度相似度的权值；

β为八连通邻域重心之间的欧氏距离的权值；

γ为八连通邻域之间的影响因子的权值；

Sim为灰度相似度；

d为八连通邻域重心之间的欧氏距离；

p为八连通邻域之间的影响因子。

(三)有益效果

本发明有如下优点：针对分水岭方法的过分割问题，引入马尔可夫随机场，利用面积和相对高度作为二维观测数据，结合图像的二维空间特征相关性，将过分割问题转化为求解图像区域标记场问题，进一步的，针对ISING模型的一维双点势能不能很好的描述邻域间的二维空间关系的问题，在将马尔可夫随机场引入分水岭变换解决过分割问题的基础上，依据邻域间的二维空间关系赋予不同势能，以抑制噪声影响，最后计算过分割区域与邻域的紧密度，选择紧密度最大的邻域并与之合并，解决过分割问题，从而实现粘连果蔬合格率的在线监测。

附图说明

图1是本发明的流程图；

图2是本发明实施例的主要步骤。

具体实施方式

下面结合附图和实施例，对本发明的具体实施方式作进一步详细描述。以下实施例用于说明本发明，但不用来限制本发明的范围。

现有对粘连果蔬分割的方法存在误判率高和耗时长等不足，有些方法甚至会对果蔬本身造成破损。因此，急需一种能够快速、准确、自动、无损地分割粘连果蔬的方法。

本发明针对现有的果蔬分割方法的不足提出了一种基于马尔科夫随机场的粘连物分割方法。包括以下步骤：

S1：将粘连物图像分割为正确分割区域和过分割区域；

S11：将粘连物图像转换为距离灰度图像；

S12：通过分水岭法对所述距离灰度图像进行处理得到初始分水岭；

S13：通过所述初始分水岭得到所述距离灰度图像的分割区域；

S14：通过所述分割区域的二维观测特征将所述分割区域划分为正确分割区域或过分割区域。

S141：获取所述分割区域的二维观测特征；

S142：根据所述二维观测特征将所述分割区域划分为标记区域；

S143：根据所述标记区域内的区域的观测特征数据计算所述标记区域内的区域的能量；

S144：重复执行步骤S141至S143，分别比较相邻两次得到的所述标记区域内的区域的能量之间的差；

S145：若所述能量之间的差小于设定阈值，则将该所述标记区域划分为正确分割区域或过分割区域。

S2：计算所述过分割区域和该所述过分割区域的八联通区域之间的紧密度，根据所述紧密度对所述过分割区域进行合并，得到粘连物的分割结果。

所述紧密度的计算公式为：

t＝αSim+βd+γp

其中：

t为紧密度；

α为灰度相似度的权值；

β为八连通邻域重心之间的欧氏距离的权值；

γ为八连通邻域之间的影响因子的权值；

Sim为灰度相似度；

d为八连通邻域重心之间的欧氏距离；

p为八连通邻域之间的影响因子。

以下通过粘连马铃薯丁分割来对本发明进行详细说明。

S1：将马铃薯丁图像分割为正确分割区域和过分割区域；

S11、将目标提取后的二值马铃薯丁图像转换为距离灰度图像；

S12、用分水岭方法浸没距离灰度图像，得到初始分水岭；

S13：通过所述初始分水岭得到所述距离灰度图像的分割区域；

S14：通过所述分割区域的二维观测特征将所述分割区域划分为正确分割区域或过分割区域；

S141：获取所述分割区域的二维观测特征；

S142：根据所述二维观测特征将所述分割区域划分为标记区域；

S143：根据所述标记区域内的区域的观测特征数据计算所述标记区域内的区域的能量；

S144：重复执行步骤S141至S143，分别比较相邻两次得到的所述标记区域内的区域的能量之间的差；

S145：若所述能量之间的差小于设定阈值，则将该所述标记区域划分为正确分割区域或过分割区域。

具体步骤如下：

令区域集记为R＝{r₁,r₂,...,r_n}，n为区域总数。令观测数据集为F＝{f₁,f₂}，其中，f₁为面积；f₂为相对高度。令为标记场，分别表示正确分割区域和过分割区域，其中，i＝1,2,...n。将过分割问题转化为求解标记场L：

$P(L=L_{r_i}|F=f)=\frac{P(L=L_{r_i})P(F=f|L=L_{r_i})}{P(F=f)},$ 其中，f＝f₁，f₂

式中：

$P(L=L_{r_i}|F=f)$ 为L的条件概率；

$P(F=f|L=L_{r_i})$ 为F的条件概率；

为L的先验概率，表示区域标记概率；

P(F)为观测数据的先验概率。

假设f中的两个特征分量：面积和相对高度相互独立，所以：

$P(L=L_{r_i}|F=f)=\frac{\underset{k=1}{\overset{2}{\mathrm{\Π}}}\left[P\left(f^k\right|L=L_{r_i})\right]P(L=L_{r_i})}{P(F=f)}$

式中，表示相对高度属于标记1的概率分布，表示面积属于标记2的概率分布。

给定F，有： $P(L=L_{r_i}|F=f)\∝\underset{k=1}{\overset{2}{\mathrm{\Π}}}\left[P\left(f^k\right|L=L_{r_i})\right]P(L=L_{r_i}).$

令E_L为E_F为E表示标记区域内的区域的能量，因为，E和logE有相同的单调性，所以对 $P(L=L_{r_i}|F=f)\∝\underset{k=1}{\overset{2}{\mathrm{\Π}}}\left[P\left(f^k\right|L=L_{r_i})\right]P(L=L_{r_i})$ 两边取对数，则有：

E＝E_L+E_F

其中：

E为标记区域内的区域的能量；

E_L为标记区域内的区域标记势能；

E_F为标记区域内的区域标记场能量。

$E_L=\underset{S}{\mathrm{\Σ}}{\mathrm{\β}}_s\×\frac{k_r}{\underset{r\∈N_s}{\mathrm{\Σ}}{k}_r}\×\underset{r\∈N_s}{\mathrm{\Σ}}\mathrm{\δ}(L_s,L_r)$

其中：

S为所有标记区域的集合；

β_s为区域s的耦合系数；

k_r为区域s与区域r之间的紧密度；

N_s为区域s的所有邻域；

δ(L_s,L_r)为势团势能V_{s,r}(L_s,L_r)的符号函数。

势团势能的计算公式为：

$V_{\{s,r\}}(L_s,L_r)=\left\{\begin{array}{c}-{\mathrm{\β}}_s\×\frac{k_r}{\underset{r\∈N_s}{\mathrm{\Σ}}{k}_r},L_s=L_r\\ {\mathrm{\β}}_s\×\frac{k_r}{\underset{r\∈N_s}{\mathrm{\Σ}}{k}_r},L_s\≠L_r\end{array}\right.$

其中：

V_{s,r}(L_s,L_r)为标记为L_s和L_r的标记场之间的势团势能；

$E_F=\underset{S,L_S=m}{\mathrm{\Σ}}\left\{\underset{k=1}{\overset{2}{\mathrm{\Σ}}}\right[\frac{{(f_s^k-{\mathrm{\μ}}_m^k)}^2}{2{\mathrm{\σ}}_m^{k2}}+\log \left(\sqrt{2\mathrm{\π}}{\mathrm{\σ}}_m^{\text{k}}\right)\left]\right\}$

其中：

S,L_S＝m为所有第m类区域的集合；

为第m类区域的第k维观测数据；

为标记区域内的区域的第k维观测数据的均值；

为标记区域内的区域的第k维观测数据的方差。

假设P(L)符合MRF模型的条件，因为对于大多数图像而言，当前区域是否为过分割区域只由其邻域决定。由Hammersley-Clifford定理：MRF与Gibbs随机场等价，因此，定义Gibbs随机场能量函数以确定MRF标记场的先验概率：

$P\left(L\right)=z^{-1}\×e^{-\frac{1}{T}{U}_1\left(L\right)}$

其中：z为归一化常数，形变系数T控制P(L)的形状，T越大，P(L)的分布就越平坦。C为双点势团的集合，V_c(L_c)为势团势能，V_c(L_c)的计算公式为：

$V_c\left(L_c\right)=V_c(L_s,L_r)=V_{\{s,r\}}(L_s,L_r)=\left\{\begin{array}{c}-{\mathrm{\β}}_s\×\frac{k_r}{\underset{r\∈N_s}{\mathrm{\Σ}}{k}_r},L_s=L_r\\ {\mathrm{\β}}_s\×\frac{k_r}{\underset{r\×N_s}{\mathrm{\Σ}}{k}_r},L_s\≠L_r\end{array}\right.$

式中：β_s为耦合系数，控制区域的同构性，在各向同性的系统中，β_s是一个常数。区域r为区域s的邻域。N_s为区域s的所有邻域，k_r为区域s与区域r的紧密度，通过计算区域r和区域s的公共边长得到，为邻域间紧密度的归一化常数。因此E_L为：

$E_L=\underset{S}{\mathrm{\Σ}}{\mathrm{\β}}_s\×\frac{k_r}{\underset{r\∈N_s}{\mathrm{\Σ}}{k}_r}\×\underset{r\∈N_s}{\mathrm{\Σ}}\mathrm{\δ}(L_s,L_r)$

其中：

S为所有标记区域的集合；

β_s为区域s的耦合系数；

k_r为区域s与区域r之间的紧密度；

N_s为区域s的所有邻域；

δ(L_s,L_r)为势团势能V_{s,r}(L_s,L_r)的符号函数。

$\mathrm{\δ}(L_s,L_r)=\left\{\begin{array}{cc}-1,& L_s=L_r\\ 1,& L_s\≠L_r\end{array}\right.$

δ(L_s,L_r)为V_{s,r}(L_s,L_r)的符号函数。

对于E_F，令和分别为第m类标号的第k个特征分量的均值和方差，表示区域s的第k个特征分量，有：

$p\left(f_s^k\right|l_s=m)=\frac{1}{\sqrt{2\mathrm{\π}}{\mathrm{\σ}}_m^k}\exp [-\frac{{(f_s^k-{\mathrm{\μ}}_m^k)}^2}{2{\mathrm{\σ}}_m^{k2}}]$

因此，E_F(f)表示为：

$E_F=\underset{S,L_S=m}{\mathrm{\Σ}}\left\{\underset{k=1}{\overset{2}{\mathrm{\Σ}}}\right[\frac{{(f_s^k-{\mathrm{\μ}}_m^k)}^2}{2{\mathrm{\σ}}_m^{k2}}+\log \left(\sqrt{2\mathrm{\π}}{\mathrm{\σ}}_m^{\text{k}}\right)\left]\right\}$

其中：

S,L_S＝m为所有第m类区域的集合；

为第m类区域的第k维观测数据；

为标记区域内的区域的第k维观测数据的均值；

为标记区域内的区域的第k维观测数据的方差。

用标记的特征向量计算第m类标号(m为标记场L中的一个标号)的均值μ_m和第m类标号的均方差σ_m。

L_S为S的标号，这里S,L_S＝m表示所有属于标号m的区域的集合S，为第m类区域的第k维观测数据。

${\mathrm{\μ}}_m^k=\frac{1}{N_m}\underset{S,L_S=m}{\mathrm{\Σ}}{f}_s^k$

${\mathrm{\σ}}_m^k=\sqrt{\frac{1}{N_m-1}\underset{S,L_S=m}{\mathrm{\Σ}}{(f_s^k-{\mathrm{\μ}}_m^k)}^2}$

其中：

N_m为第m类区域的特征向量个数；

为第m类区域的第k维观测数据；

S为所有标记区域的集合；

S,L_S＝m为所有第m类区域的集合；

L_S为S的标号。

用模拟退火法最小化下式并重新标记图像，更新μ_m和σ_m，直到相邻两次的差值小于设定阈值。得到标记场L；

L＝argmin_l(E_L+E_F)

通过标记场对分割区域进行正确分割区域和过分割区域标识。

S2：计算马铃薯丁的过分割区域和过分割区域的八联通区域之间的紧密度，根据所述紧密度对所述过分割区域进行合并，得到马铃薯丁的的分割结果。

原始图像中八连通邻域的灰度直方图和面积比例的相似度计算灰度相似度，记为Sim：

$\mathrm{Sim}=\frac{\underset{j=1}{\overset{256}{\mathrm{\Σ}}}j\×N_j}{\underset{j=1}{\overset{256}{\mathrm{\Σ}}}{N}_j}/\frac{\underset{i=1}{\overset{256}{\mathrm{\Σ}}}i\×N_i}{\underset{i=1}{\overset{256}{\mathrm{\Σ}}}{N}_i}$

其中：

N_j为灰度值j对应的像素个数；

N_i为灰度值i对应像素个数；

j为当前过分割区域的八连通邻域的灰度值；

i为当前过分割区域的灰度值。

距离为八连通邻域重心之间的欧氏距离，记为d；八连通邻域之间的影响因子通过归一化公共边长获得，记为p：

$p=\frac{{\mathrm{Pub}}_r}{\underset{r\∈N_s}{\mathrm{\Σ}}{k}_r}$

式中，Pub_r为八连通邻域的公共边长，N_s为当前区域的所有邻域。

因此，过分割区域与八连通邻域的紧密度t为：

t＝αSim+βd+γp

其中：

t为紧密度；

α为灰度相似度的权值；

β为八连通邻域重心之间的欧氏距离的权值；

γ为八连通邻域之间的影响因子的权值；

Sim为灰度相似度；

d为八连通邻域重心之间的欧氏距离；

p为八连通邻域之间的影响因子。

在流水线上采集一帧粘连马铃薯丁图像，如图2，

本发明对瓜果、块状农产品的分级具有通用性，但由于产品种类很多，因此本发明只举一个用于粘连鲜切马铃薯丁分割的实施实例，其它瓜果和块状农产品的切块粘连分割也可以参照该实施实例的方法，具体针对所测瓜果、块状农产品的切块特征，改变相关参数，就可以对新的切块粘连图像进行分割。

以上实施方式仅用于说明本发明，而并非对本发明的限制，有关技术领域的普通技术人员，在不脱离本发明的精神和范围的情况下，还可以做出各种变化和变型，因此所有等同的技术方案也属于本发明的范畴，本发明的专利保护范围应由权利要求限定。

Claims (7)

1.一种基于马尔科夫随机场的粘连物分割方法，其特征是，该方法包括以下步骤：

S1：将粘连物图像分割为正确分割区域和过分割区域；

S2：计算所述过分割区域和该所述过分割区域的八连通区域之间的紧密度，根据所述紧密度对所述过分割区域进行合并，得到粘连物的分割结果；

所述步骤S1具体为：

S11：将粘连物图像转换为距离灰度图像；

S12：对所述距离灰度图像进行处理得到初始分水岭；

S13：通过所述初始分水岭得到所述距离灰度图像的分割区域；

S14：通过所述分割区域的二维观测特征将所述分割区域划分为正确分割区域或过分割区域；

所述步骤S14具体为：

S141：获取所述分割区域的二维观测特征；

S142：根据所述二维观测特征将所述分割区域划分为标记区域；

S143：根据所述标记区域内的区域的观测特征数据计算所述标记区域内的区域的能量；

S144：重复执行步骤S141至S143，分别比较相邻两次得到的所述标记区域内的区域的能量之间的差；

S145：若所述能量之间的差小于设定阈值，则将该所述标记区域划分为正确分割区域或过分割区域。

2.如权利要求1所述的方法，其特征是，所述步骤S12具体为：通过分水岭法对所述距离灰度图像进行处理得到初始分水岭。

3.如权利要求1所述的方法，其特征是，所述标记区域内的区域的能量的计算公式为：

E＝E_L+E_F

其中：

E为标记区域内的区域的能量；

E_L为标记区域内的区域标记势能；

E_F为标记区域内的区域标记场能量。

4.如权利要求3所述的方法，其特征是，所述的E_L计算公式为：

$E_L=\underset{S}{\mathrm{\Σ}}{\mathrm{\β}}_s\×\frac{k_r}{\underset{r\∈N_s}{\mathrm{\Σ}}{k}_r}\×\underset{r\∈N_s}{\mathrm{\Σ}}\mathrm{\δ}(L_s,L_r)$

其中：

S为所有标记区域的集合；

β_s为区域s的耦合系数；

k_r为区域s与区域r之间的紧密度；

N_s为区域s的所有邻域；

δ(L_s,L_r)为势团势能V_{s,r}(L_s,L_r)的符号函数。

5.如权利要求4所述的方法，其特征是，所述势团势能的计算公式为：

$V_{\{s,r\}}(L_s,L_r)=\left\{\begin{array}{c}-{\mathrm{\β}}_s\×\frac{k_r}{\underset{r\∈N_s}{\mathrm{\Σ}}{k}_r},L_s=L_r\\ {\mathrm{\β}}_s\×\frac{k_r}{\underset{r\∈N_s}{\mathrm{\Σ}}{k}_r},L_s\≠L_r\end{array}\right.$

其中：

V_{s,r}(L_s,L_r)为标记为L_s和L_r的标记场之间的势团势能。

6.如权利要求3所述的方法，其特征是，所述的E_F计算公式为：

$E_F=\underset{S,L_S=m}{\mathrm{\Σ}}\left\{\underset{k=1}{\overset{2}{\mathrm{\Σ}}}\right[\frac{{(f_s^k-{\mathrm{\μ}}_m^k)}^2}{2{\mathrm{\σ}}_m^{k2}}+\log \left(\sqrt{2\mathrm{\π}}{\mathrm{\σ}}_m^k\right)\left]\right\}$

其中：

S,L_S＝m为所有第m类区域的集合；

为第m类区域的第k维观测数据；

为标记区域内的区域的第k维观测数据的均值；

为标记区域内的区域的第k维观测数据的方差。

7.如权利要求1所述的方法，其特征是，所述紧密度的计算公式为：

t＝αSim+βd+γp

其中：

t为紧密度；

α为灰度相似度的权值；

β为八连通邻域重心之间的欧氏距离的权值；

γ为八连通邻域之间的影响因子的权值；

Sim为灰度相似度；

d为八连通邻域重心之间的欧氏距离；

p为八连通邻域之间的影响因子。