[发明专利]一种基于多尺度水平集的眼底血管图像分割系统及方法

权利要求书

1.一种基于多尺度水平集的眼底血管图像分割系统，其特征在于：包括输入模块、预处理模块、尺度构建模块、水平集模块以及输出模块；

所述输入模块，用于接收用户输入的眼底血管图像，然后将其传递到预处理模块；

所述预处理模块，用于对接收的眼底血管图像进行扩展，扩展的结果传递到尺度构建模块；

所述尺度构建模块，用于对眼底血管图像找到最佳的血管响应尺度，进而基于该尺度构建高斯曲线，将该高斯曲线与图像卷积的结果减去原图像得到图像M，将图像M传递到水平集模块；

所述水平集模块用于对图像M采用水平集算法进行血管分割，并将分割结果传递到输出模块；

所述输出模块用于对分割结果图像进行输出显示。

2.一种基于多尺度水平集的眼底血管图像分割方法，通过权利要求1所述一种基于多尺度水平集的眼底血管图像分割系统实现，其特征在于：包括以下步骤：

步骤1：输入模块接收待分割眼底血管图像A，图像A的宽度为W，高度为H；

步骤2：预处理模块对接收的待分割眼底血管图像A进行构造，生成无边缘背景的眼底血管图像I，其分辨率与图像A的分辨率相同；将图像A转为灰度图像J，图像J的全部区域设为Ω，若待分割区域充满Ω，则直接获得图像P，若待分割区域没有充满Ω，对图像内的有效边界区域进行镜像扩展使整个图像内的黑色填充区被覆盖，扩展后的图像为P，同时，I＝P，执行步骤2.1；否则I＝J，执行步骤3；

步骤2.1：确定阈值Q将灰度图像J划分为前景R和背景B，其中R和B是二值图像，前景区域处R＝1，B＝0，背景区域处R＝0，B＝1；

步骤2.2：构造前景区域边界图像R₂，将前景R进行腐蚀操作得到R₁，R₂＝R-R₁，其中图像R₂为二值图像；

步骤2.3：初始化变量p＝1，其中P是外循环中像素点的初始化横坐标值；

步骤2.4：如果p≤W，且外循环中像素点的初始化纵坐标值q＝1时,执行步骤2.5，否则，执行步骤3；

步骤2.5：如果q≤H，执行步骤2.6，否则,令p＝p+1,执行步骤2.4；

步骤2.6：如果当前点(p,q)满足B＝1，执行步骤2.7，否则令q＝q+1，执行步骤2.5；

步骤2.7：初始化变量m＝1，m为图像内循环中像素点的初始横坐标值，初始化其中r为初始化图像中任意两个像素之间的距离，初始坐标(a，b)为两个随机正整数；

步骤2.8：如果m≤W，变量n＝1,n为内循环中像素点的纵坐标的初始值，执行步骤2.9，否则，执行步骤2.13；

步骤2.9：如果n≤H，执行步骤2.10，否则,令m＝m+1,执行步骤2.8；

步骤2.10：如果点(m,n)满足R₂＝1，执行步骤2.11，否则,令n＝n+1,执行步骤2.9；

步骤2.11：计算点(p,q)与点(m,n)之间的距离dist((p,q),(m,n))，其中dist为距离度量函数，若r＜dist((p,q),(m,n))，则r＝dist((p,q),(m,n))，同时令a＝m，b＝n；

步骤2.12：令n＝n+1，跳转到步骤2.9；

步骤2.13：点(p,q)以点(a,b)为中心，r为对称距离，得到R中一个对应像素点L，将像素点L的值赋值给当前点Z；其中L_x＝2a-p，L_y＝2b-q，其中L_x是像素点L的横坐标L_y是像素点L的纵坐标；

步骤2.14：令q＝q+1,跳转到步骤2.5；

步骤3：对参数最小尺度σ_min、最大尺度σ_max、常数系数γ、常数系数c、尺度步长Δσ进行设定，其中σ_min＜σ_max；

步骤4：对初始尺度σ_i进行设定，使σ_i＝σ_min；

步骤5：构造图像I的Hessian矩阵，通过最佳响应尺度σ构建高斯函数，将图像I与高斯函数进行卷积操作的结果减去图像I得到图像M，具体步骤如下：

步骤5.1：初始化变量p＝1，尺度变量κ＝1，变量ν_max＝0，其中v_max是像素的最大响应值变量；

步骤5.2：如果p≤W，则变量q＝1，执行步骤5.3，否则，执行步骤6；

步骤5.3：如果q≤H，执行步骤5.4，否则，令p＝p+1，执行步骤5.2；

步骤5.4：判断当前尺度σ_i是否大于最大尺度σ_max，若大于σ_max，则令σ_i＝σ_min，执行步骤5.12，否则，执行步骤5.5；

步骤5.5：构造点(p,q)的Hessian矩阵：

其中，

D＝∫∫I(p-u,q-v)×A(u,v)dudv

E＝∫∫I(p-u,q-v)×B(u,v)dudv

F＝∫∫I(p-u,q-v)×C(u,v)dudv

A(u,v)、B(u,v)和C(u,v)是基于尺度σ_i构造高斯二阶偏导数

其中D、E、F为图像I与对应的高斯二阶偏导数A、B、C的卷积；

步骤5.6：计算步骤5.5中得到的Hessian矩阵的两个特征值λ₁、λ₂，其中

其中

步骤5.7：计算在当前尺度σ_i下，像素(p,q)的响应值

其中

式中，R_B是线性形状的各向异性程度，S是范数；

步骤5.8：若则令且κ＝σ_i；

步骤5.9：基于尺度κ构建一个高斯核得到高斯曲线g_κ；

步骤5.10：将步骤5.9得到的高斯曲线g_κ与图像I进行卷积操作，将卷积结果减去图像I得到图像M；

M(p,q)＝(g_κ*I)(p,q)-I(p,q)

步骤5.11：更新尺度σ_i＝σ_i+Δσ，Δσ为是尺度步长，跳转到步骤5.4；

步骤5.12：更新坐标q＝q+1，跳转到步骤5.3；

步骤6:在图像M上设定初始的零水平集分割曲线，根据设定的零水平集生成高一维的水平集函数φ⁰；

所述设定的初始零水平集分割曲线是指定形状的闭合曲线是随机生成的闭合曲线，所述水平集函数φ是根据设定的初始水平集分割曲线生成的高一维符号距离函数；

步骤7:建立水平集分割的能量泛函E(φ)，并设定水平集方法中的加权系数α₁、加权系数α₂、曲线平滑项系数v、面积平滑项系数μ；

所述能量泛函如下所示：

E(φ)＝E_D+vL(φ)+μP(φ)

其中，E_D为图像数据项，L(φ)为长度约束项，P(φ)为面积约束项，H(φ(p,q))为φ(p,q)的Heaviside函数；

E_D(φ)＝α₁∫_Ω|M(p,q)-c₁|²H(φ(p,q))dpdq+α₂∫_Ω|M(p,q)-c₂|²(1-H(φ(p,q)))dpdq

P(φ)＝∫H(φ(p,q))dpdq

其中c₁是零水平集分割曲线内的灰度均值，c₂是零水平集分割曲线外的灰度均值；

步骤8：设定最大迭代次数J和图像水平集的能量泛函E(φ)的收敛阈值ε；

步骤9：分别更新图像水平集函数集合H(φ^j)、1-H(φ^j)、C^j、φ^j、得到第j次迭代的能量泛函E^j；

具体地：

其中Δt为时间步长，为迭代更新的水平集函数的时间偏导数，所述时间偏导数公式如下所示：

步骤10：对是否达到迭代终止条件进行判断：若当前迭代次数j达到最大迭代次数J或者能量泛函E(φ)相邻两次的迭代结果差值小于提前设定的收敛阈值ε，则执行步骤11，否则，迭代次数加1，返回步骤9；

步骤11：根据当前更新后的图像水平集函数φ^j构造图像，即图像的分割结果。