技术领域
[0001] 本
发明属于图像分割技术领域，尤其涉及一种基于中心点检测的CBCT牙齿分割方法。
背景技术
[0002] 随着社会发展，越来越多的成年人在寻求牙齿正畸
治疗，以改善他们的笑容，改正牙齿咬合状况或纠正因受伤、
疾病或长久忽视
口腔护理所造成的其它问题。为全面、精确、详实地获取患者数据，进而为口腔疾病的诊断和治疗提供有
力保障，医学成像技术成为不可或缺的技术手段。口腔CT设备是牙科的革命性工具，CBCT自上世纪九十年代后期应用于口腔医学领域，因其能够真实反映颌面部的三维结构，克服二维成像的弊端，可以更加立体的评估牙齿和骨骼的关系，使正畸医生能够制定更合理的方案。但由于CBCT牙齿切片图像
亮度低、
对比度差、噪点多、边界不明显、结构异常复杂，个体呈多样性；并且CT因受噪声、强度非均匀性产生的伪影及不同组织强度相同的影响使得牙齿的自动分割非常困难，且大多数图像灰度不均匀。如何处理有效的分割牙齿图像成为了获取牙齿信息的关键步骤。
[0003] 目前主要的牙齿分割
算法有：手动分割，其缺点是需要多人参与、分割时间长、测量
精度受到操作者操作影响；自适应
阀值分割，但因为存在噪点多、相邻牙齿边界模糊和牙齿灰度值不均匀这些问题，极易受到非有效信息的干扰；形态学分割，由于对边缘过于敏感，经常出现过度分割；基于活动轮廓模型的
水平集分割，如DRLSE、MICO存在的问题是对于初始化要求过高需要依赖于较为精准先验知识，而且缺少对演化表面的局部控制。由于医学图像的成像技术
缺陷和牙齿结构的特殊性，牙齿的分割处理目前都要借助一定量的交互操作才能实现。因此，引入形状先验信息来简化分割问题是必不可少的。为了实现ROI的精确自动分割，本文考虑利用检测到的中心点信息作为分割算法的先验信息，提出基于中心点检测的牙齿图像分割算法，使分割更准确，鲁棒性更好。基于中心点检测的牙齿图像分割算法充分利用检测到的中心点信息，进行初始层的粗分割，并利用粗分割结果代替手动初始化，通过关键区域的先验信息直接融入到待分割图像中，实现一种全自动的“专家”指导下的图像分割。
发明内容
[0004] 针对
现有技术中的上述不足，本发明提供的一种基于中心点检测的CBCT牙齿分割方法解决了针对现有的大部分牙齿分割算法都需要手动初始化，结合提高牙齿分割精度的问题。
[0005] 为了达到以上目的，本发明采用的技术方案为：
[0006] 本方案提供一种基于中心点检测的CBCT牙齿图像分割方法，包括以下步骤：
[0007] S1、预处理：获取牙齿的原始图像，并根据所述牙齿的原始图像利用MIP算法计算得到包裹牙齿的区域大小；
[0008] S2、粗分割：根据所述包裹牙齿的区域大小选取合适的初始层，并利用分水岭算法计算得到牙齿所述初始层的粗分割结果；
[0009] S3、双水平集精分割：将所述粗分割结果作为初始化，并利用双水平集DRLSE模型对牙齿初始层进行精分割处理，得到上下两排牙齿初始层的二维分割结果；
[0010] S4、层间
迭代分割：利用层间信息将所述上下两排牙齿初始层的二维分割结果进行最优
阈值处理，并根据其处理结果利用双水平集DRLSE模型进行向上或向下逐层迭代得到CBCT图像牙齿每一层的二维分割结果；
[0011] S5、输出三维牙齿结构：利用DRLSE模型对所述CBCT图像牙齿的二维分割结果进行分割处理，得到CBCT牙齿图像的三维分割结果，从而完成CBCT牙齿图像的分割。
[0012] 进一步地，所述步骤S1包括以下步骤：
[0013] S101、获取牙齿原始图像，以DICOM格式保存，并逐层读取牙齿切片图像；
[0014] S102、对所述牙齿原始图像进行分段线性变换处理，将牙齿原始图像的灰度归一化到[0,255]；
[0015] S103、利用MIP投影算法对经灰度归一化处理后的图像在x，y和z三个方向分别进行投影，得到包裹牙齿的区域大小；
[0016] 所述对图像的x，y和z三个方向分别进行投影的表达式如下：
[0017] xmip(j,k)＝max(xmip(j,k),a(i,j,k))；
[0018] ymip(i,k)＝max(ymip(i,k),a(i,j,k))；
[0019] zmip(i,j)＝max(zmip(i,j),a(i,j,k))；
[0020] 其中，xmip(j,k)表示图像x方向的投影，ymip(i,k)表示图像y方向的投影，zmip(i,j)表示图像z方向的投影，i表示范围从1到x方向的图像大小，j表示范围从1到y方向的图像大小，k表示范围从1到z方向的图像大小，a(i,j,k)表示图像在i，j和k处的灰度大小。
[0021] 再进一步地，所述步骤S2包括以下步骤：
[0022] S201、根据所述包裹牙齿的区域大小分别选取上下两排牙齿的初始层；
[0023] S202、根据所述上下两排牙齿的初始层分别检测上下两排牙齿的中心点，并将所述中心点作为前景，预设的阈值作为背景；
[0024] S203、利用分水岭算法将所述前景和背景分别进行内外部标记，分别得到上下两排牙齿的粗分割结果。
[0025] 再进一步地，所述步骤S3包括以下步骤：
[0026] S301、根据所述步骤S2得到的粗分割结果，按照牙齿的交替顺序对上下两排牙齿进行双水平集的两个初始化处理，得到上下两排牙齿初始层的初始化；
[0027] S302、利用双水平集DRLSE模型分别对上下两排牙齿初始层进行精分割，从而得到上下排牙齿初始层的二维分割结果。
[0028] 再进一步地，所述步骤S4包括以下步骤：
[0029] S401、利用层间信息将所述上下排牙齿初始层的二维分割结果进行最优阈值处理，并将处理后的结果作为向上一层或向下一层的初始化；
[0030] S402、利用双水平集DRLSE模型分别对上下两排牙齿依次进行迭代分割处理，从而得到CBCT图像牙齿每一层的二维分割结果。
[0031] 再进一步地，所述步骤S401中进行最优阈值处理的表达式如下：
[0032]
[0033] 其中，xi表示[0,255]的
像素灰度，yi表示该灰度在图像曲线内部的像素点个数，A表示幅值，S表示曲线宽度，e表示自然常数，μ表示控制演化的参数。
[0034] 再进一步地，所述步骤S402具体为：
[0035] 利用双水平集DRLSE模型分别对上下两排牙齿先从初始层到
牙冠方向进行逐层分割，再从初始层到
牙根方向逐层分割，从而得到CBCT图像牙齿的二维分割结果。
[0036] 再进一步地，所述步骤S302和步骤S402中每次迭代分割过程中的两条双水平集演化曲线需满足如下条件：
[0037] φ1＝max(φ1,-φ2)
[0038] φ2＝max(-φ1,φ2)
[0039] 其中，φ1表示第一条水平曲线，φ2表示第二条水平曲线，max(·)表示取最大值。
[0040] 再进一步地，所述步骤S302和步骤S402利用双水平集DRLSE模型分割的表达式如下：
[0041]
[0042]
[0043]
[0044]
[0045]
[0046] 其中， 表示均表示时间演化方程，μ、λ以及α均表示控制演化的参数，δ表示单位冲激函数，div(·)表示散度，φ1表示第一条水平曲线，φ2表示第二条水平曲线，表示两条水平集对应的梯度，dp(s)表示参数s定义的函数，p'(s)表示对函数p的求导，s表示[0,1]的空间参数，p表示距离正则化项定义的函数，p2(·)表示关于水平集梯度模的函数， 表示水平集曲线的梯度模， 表示梯度算子，Gσ表示标准差为σ的高斯函数，I表示整幅图像，*表示卷积算子，g表示
边缘检测因子。
[0047] 再进一步地，所述步骤S5具体为：
[0048] S501、将上下两排牙齿的
单层CBCT牙齿图像的二维分割结果进行堆叠，并将堆叠结果作为三维DRLSE模型的初始化；
[0049] S502、利用三维DRLSE模型对所述堆叠结果进行迭代处理，得到CBCT牙齿图像的三维分割结果，从而完成CBCT牙齿图像的分割。
[0050] 本发明的有益效果：
[0051] (1)本发明将检测到的牙齿中心点当作粗分割算法的先验信息，利用前景区域的先验信息得到图像里面牙齿的粗略边界，并将其作为水平集模型的初始化。借此可以降低医生手工干预造成的时间浪费和人工初始化不合理问题的情况发生；
[0052] (2)本发明借助层与层之间的信息进行约束的双水平集函数的操作能够有效降低函数自身的噪声敏感性，防止相邻牙齿吸引造成的分割不准确的情况，并可充分利用牙齿区域信息，将全口腔牙齿进行双水平集进行有效率的分割，得到更好的分割结果。
附图说明
[0054] 图2为本
实施例中的牙齿CBCT影像图。
[0055] 图3为本实施例中粗分割结果示意图
[0056] 图4为本实施例中双水平集单层分割的结果图。
[0057] 图5为本实施例中牙齿分割过程的示意图。
[0058] 图6为本实施例中曲线内部灰度高斯拟合示意图。
具体实施方式
[0059] 下面对本发明的具体实施方式进行描述，以便于
本技术领域的技术人员理解本发明，但应该清楚，本发明不限于具体实施方式的范围，对本技术领域的普通技术人员来讲，只要各种变化在所附的
权利要求限定和确定的本发明的精神和范围内，这些变化是显而易见的，一切利用本发明构思的发明创造均在保护之列。
[0060] 实施例
[0061] 本发明针对现有的大部分牙齿分割算法都需要手动初始化，结合提高牙齿分割精度的情况。本发明提出一种有效提高自动化程度的思路，并在整个输出三维牙齿结构分割的过程中引入双水平集及约束处理，便于更快、更方便也更准确的得到牙齿三维结果。本文考虑利用检测到的中心点信息作为分割算法的先验信息，提出基于中心点检测的牙齿图像分割算法，使分割更准确，鲁棒性更好。基于中心点检测的牙齿图像分割算法充分利用检测到的中心点信息，进行初始层的粗分割，并利用粗分割结果代替手动初始化，通过引入关键区的先验信息直接融入到待分割图像中，实现一种全自动的“专家”指导下的图像分割。本发明针对现有的大部分牙齿分割算法都需要手动初始化，结合提高牙齿分割精度的情况。
[0062] 如图1所示，本发明公开了一种基于中心点检测的CBCT牙齿图像分割方法，其实现方法如下：
[0063] S1、预处理：获取牙齿原始图像，并根据所述牙齿原始图像利用MIP算法计算得到包裹牙齿的区域大小，其实现方法如下：
[0064] S101、获取牙齿原始图像，以DICOM格式保存，并逐层读取牙齿切片图像；
[0065] S102、对所述牙齿原始图像进行分段线性变换处理，将牙齿原始图像的灰度归一化到[0,255]；
[0066] S103、利用MIP投影算法对经灰度归一化处理后的图像在x，y和z三个方向分别进行投影，得到包裹牙齿的区域大小。
[0067] 本实施例中，首先获取原始图像，以标准DICOM格式保存，逐层读取牙齿切片图像。并求取原始图像的灰度直方图，进行分段线性变换将图像灰度归一化到[0,255]，然后利用MIP算法将原始
图像处理得到包裹牙齿区域大小，减小处理图像大小，便于提高分割速度。
[0068] 本实施例中，如图2所示，所述步骤S1中对得到原始图像的灰度直方图显示图像灰度范围较大近3000，为了便于图像处理，进行灰度变换到[0,255]，然后由于单层图片太大会造成处理速度变慢，因此，利用MIP投影算法对x，y和z三个方向分别进行投影，原始图像大小eg:512×512×512处理得到包裹牙齿区域大小eg:250×300×280)即可。最大
密度投影(Maximum Intensity Projection,MIP)：顾名思义，把最大的值投射出来。即，计算每一条光线穿过的像素
块中灰度最大的，作为输出结果。对应x，y和z三个方向情况：
[0069] xmip(j,k)＝max(xmip(j,k),a(i,j,k))；
[0070] ymip(i,k)＝max(ymip(i,k),a(i,j,k))；
[0071] zmip(i,j)＝max(zmip(i,j),a(i,j,k))；
[0072] 其中，xmip(j,k)表示图像x方向的投影，ymip(i,k)表示图像y方向的投影，zmip(i,j)表示图像z方向的投影，i表示范围从1到x方向的图像大小，j表示范围从1到y方向的图像大小，k表示范围从1到z方向的图像大小，a(i,j,k)表示图像在i，j和k处的灰度大小。
[0073] S2、粗分割：根据所述包裹牙齿的区域大小选取合适的初始层，并利用分水岭算法计算得到牙齿所述初始层的粗分割结果，其实现方法如下：
[0074] S201、根据所述包裹牙齿的区域大小分别选取上下两排牙齿的初始层；
[0075] S202、根据所述上下两排牙齿的初始层分别检测上下两排牙齿的中心点，并将所述中心点作为前景，预设的阈值作为背景；
[0076] S203、利用分水岭算法将所述前景和背景分别进行内外部标记，分别得到上下两排牙齿的粗分割结果。
[0077] 本实施例中，首先选取合适的初始层；利用初始层检测得到的中心点当作前景，阈值作为背景，利用控制标记分水岭算法得到初始层的粗分割结果，确保可以得到一个比较好的牙齿不黏连的初始化。
[0078] 本实施例中，所述步骤S2中牙齿之间会有较多的粘连，为了区分开牙齿粘连的可能性，我们尽可能选取牙冠层到有
牙周组织出现层的
中间层作为初始层。但由于要全部考虑到人体单层16颗牙齿的情况，初始层确保每颗牙齿尽可能互相不黏连；因此我们在进行精分割之前，首先利用初始层检测得到的中心点当作前景，阈值作为背景，利用控制标记分水岭算法得到初始层的粗分割结果，确保可以得到一个比较好的牙齿不黏连的初始化。分水岭变换得到的是输入图像的集水盆图像，集水盆之间的边界点，即为分水岭。显然，分水岭表示的是输入图像极大值点。因此，为得到图像的边缘信息，通常把梯度图像作为输入图像，即：
[0079] g(x,y)＝grad(f(x,y))
[0080] 其中，f(x,y)表示原始图像，grad(·)表示梯度运算。但由于噪声和梯度的其他不规则性造成的过度分割，因此引入标记，将前景和背景分别进行内外部标记。如图3所示，即可得到较好的牙齿粗分割结果。
[0081] S3、双水平集精分割：如图4所示，将所述粗分割结果作为初始化，并利用双水平集DRLSE模型对牙齿初始层进行精分割处理，得到上下两排牙齿初始层的二维分割结果，其实现方法如下：
[0082] S301、根据所述步骤S2得到的粗分割结果，按照牙齿的交替顺序对上下两排牙齿进行双水平集的两个初始化处理，得到上下两排牙齿初始层的初始化；
[0083] S302、利用双水平集DRLSE模型分别对上下两排牙齿初始层进行精分割，从而得到上下排牙齿初始层的二维分割结果。
[0084] 本实施例中，利用DRLSE模型得到初始层的牙齿分割结果，并利用层间信息，把上一层的结果当作下一层的初始化，依次进行迭代分割。由于对于磨牙
位置，单个水平集方法不能很好的解决同一颗牙齿的分根问题，因此，在层与层之间的迭代分割过程中对上一层的分割结果进行最优阈值方法进行优化处理。另外利用单个水平集对每颗牙齿分割很可能会造成两颗相近牙齿分割结果被旁边牙齿吸引过去的情况，为了保证单颗牙齿的分割精度，我们拟采用双水平集方法进行分割处理，在每一次迭代过程中对双水平集进行交并处理，保证两条曲线在演化过程中不会出现相互交错或者融合成一条边界的情况，有利于分割结果的提高。通过迭代的方式求解水平集函数即可得到上下两排牙每一层分割结果。
[0085] 本实施例中，所述步骤S3主要思路是当对上下两排牙齿分别选取合适的初始层，下排牙从初始层到牙根分割，从初始层到牙冠分割；上排牙从初始层到牙冠分割，再从初始层到牙根分割，分割过程示意图如图5所示得到粗分割结果，最终将各个层的分割结果进行
三维重建。
[0086] 精分割模型利用双水平集DRLSE模型，双水平集DRLSE模型如下：
[0087]
[0088]
[0089]
[0090]
[0091]
[0092] 其中， 表示均表示时间演化方程，μ、λ以及α均表示控制演化的参数，δ表示单位冲激函数，div(·)表示散度，φ1表示第一条水平曲线，φ2表示第二条水平曲线，表示两条水平集对应的梯度，dp(s)表示参数s定义的函数，p'(s)表示对函数p的求导，s表示[0,1]的空间参数，p表示距离正则化项定义的函数，p2(·)表示关于水平集梯度模的函数， 表示水平集曲线的梯度模， 表示梯度算子，Gσ表示标准差为σ的高斯函数，I表示整幅图像，*表示卷积算子，g表示边缘检测因子。
[0093] 本实施例中，两条双水平集演化曲线需保证：
[0094] φ1＝max(φ1,-φ2)
[0095] φ2＝max(-φ1,φ2)
[0096] 其中，φ1表示第一条水平曲线，φ2表示第二条水平曲线，max(·)表示取最大值，即可保证在水平集函数迭代演化过程中两条水平集函数相互排斥，从而达到分离开两颗相邻牙齿的目的，则全口腔每一层的16颗牙齿即可按照单双数交替的形式完成两条水平集的分割，得到单层CBCT牙齿图像的二维分割结果。
[0097] 层初始化阈值处理：如图6所示，将曲线内部灰度进行高斯拟合，将数据点拟合为高斯函数类型，即使如下公式：
[0098]
[0099] 来对曲线内的像素点数据(xi,yi),(i＝1,2,3,...)，其中，xi表示[0,255]的像素灰度，yi表示该灰度在图像曲线内部的像素点个数，A表示幅值，S表示曲线宽度，e表示自然常数，μ表示控制演化的参数。
[0100] 本实施例中，对前一层精分割结果曲线内部灰度分布进行高斯拟合，根据3σ准则选取T＝μ-3σ作为当前层结果内部阈值处理，提高分割精度，减少周边区域的干扰，或者因水平集迭代出现被邻牙吸引过去的情况。
[0101] S4、层间迭代分割：利用层间信息将所述上下两排牙齿初始层的二维分割结果进行最优阈值处理，并根据其处理结果利用双水平集DRLSE模型进行向上或向下逐层迭代得到CBCT图像牙齿每一层的二维分割结果，其实现方法如下：
[0102] S401、利用层间信息将所述上下排牙齿初始层的二维分割结果进行最优阈值处理，并将处理后的结果作为向上一层或向下一层的初始化；
[0103] S402、利用双水平集DRLSE模型分别对上下两排牙齿依次进行迭代分割处理，从而得到CBCT图像牙齿每一层的二维分割结果；
[0104] 所述步骤S401中进行最优阈值处理的表达式如下：
[0105]
[0106] 其中，xi表示[0,255]的像素灰度，yi表示该灰度在图像曲线内部的像素点个数，A表示幅值，S表示曲线宽度，e表示自然常数，μ表示控制演化的参数，然后根据3σ准则进行阈值选取，其中σ为标准差；
[0107] 所述步骤S402具体为：
[0108] 利用双水平集DRLSE模型分别对上下两排牙齿先从初始层到牙冠方向进行逐层分割，再从初始层到牙根方向逐层分割，从而得到CBCT图像牙齿的二维分割结果。
[0109] S5、输出三维牙齿结构：利用DRLSE模型对所述CBCT图像牙齿的二维分割结果进行分割处理，得到CBCT牙齿图像的三维分割结果，从而完成CBCT牙齿图像的分割，其实现方法如下：
[0110] S501、将上下两排牙齿的单层CBCT牙齿图像的二维分割结果进行堆叠，并将堆叠结果作为三维DRLSE模型的初始化；
[0111] S502、利用三维DRLSE模型对所述堆叠结果进行迭代处理，得到CBCT牙齿图像的三维分割结果，从而完成CBCT牙齿图像的分割。
[0112] 本实施例中，根据得到的单层分割结果作为三维DRLSE模型的初始化，通过迭代的方式求解三维水平集函数以得到准确的三维分割结果。
[0113] 本实施例中，所述步骤S5将步骤S4中得到的二维结果堆叠起来作为三维DRLSE的初始化，利用三维DRLSE模型进行最终的三维精确分割处理，可以优化分割结果，填补二维堆叠结果带来的小孔洞、小间隙，使最终的牙齿三维分割结果更精确。
