CN104732516A - 一种基于随机方向直方图比的双阈值的血管图像处理方法 - Google Patents

Abstract

本发明公开了一种基于随机方向直方图比的双阈值的血管图像处理方法，属于图像处理技术领域，提出了一种基于随机探针检测血管图像子窗口内的图像性质的指标——方向直方图比，根据方向直方图比的值大小，能够鉴别血管图像中的噪声污染区域与血管区域，对噪声污染区域采用高阈值进行阈值化处理，从而能够尽量的消除噪声像素，而对血管区域则采用低阈值进行阈值化处理，从而尽量的保留血管像素，其中，高、低两个阈值则通过三级Otsu算法获得，通过该方法获得的血管检测结果，能够有效得去除血管图像中的噪声。

Description

一种基于随机方向直方图比的双阈值的血管图像处理方法

技术领域

本发明属于图像处理技术领域，具体涉及一种基于随机方向直方图比的双阈值的血管图像处理方法。

背景技术

为了能够理解复杂的血管图像数据，许多现有研究已经提出了各种形态学方法对血管进行量化，包括血管长度分布、血管半/直径分布、血管方向、血管区域或空间密度、血管分枝节点密度、血管分枝角度、血管端点密度和分形维数等。显然，这些度量的精确程度严重依赖于图像的血管检测输出结果。然而，由于背景着色、血管漏洞或者通过血管壁渗透到周围组织的荧光染料的存在，血管造影图像数据往往会出现部分区域的噪声污染。在这种情况下，从血管图像中进行有效的血管检测具有较大的难度。

现有研究针对血管图像数据处理问题已经提出了多种方法，根据血管图像数据所采用信息和特征种类的不同，血管图像数据处理算法主要可以三类，包括基于颜色特征(亮度)的算法，基于形状(几何特征)的算法，以及联合使用了几何特征以及亮度特征的算法。基于颜色特征的方法主要通过研究颜色聚类、直方图分割、图像熵、像素局部邻域内的统计特征、或者图像的某些属性(如边缘特征等)，进行血管的分割提取。其中，最为直接应用也最为广泛的一类算法是基于亮度信息的阈值化处理。阈值化处理根据不同的分类标准，可以进行不同的归类。根据阈值化处理时所参考信息的范围，可以分为全局方法和局部方法两类。此外，根据所采用信息种类的不同，阈值化算法又可以分为基于直方图形状的算法，基于聚类的算法，基于熵的算法，基于属性相似性的算法，基于空间分布的算法，以及基于局部统计特征的算法。基于直方图的方法根据灰度直方图的形状，特别是峰与谷分割来计算阈值。基于聚类的方法， 根据像素间的相似性或者差异性来区分前景与背景。基于熵的方法，通过最大化二值图像的熵，或者最小化原始图像与所获得的二值图像之间的交叉熵，来确定阈值。基于属性的方法，通过比较原始图像与结果二值图像之间的特征差异来确定最佳二值图像，其所比较的特征包括边缘、形状、灰度矩等。这一类算法试图在最大程度上保留原始图像的几何特征。基于空间分布的算法，主要根据相邻像素之间的相关性来确定阈值，包括共生矩阵以及二阶矩等统计相关特征。其主要的策略是，将结果二值图像与原始图像之间的分布变化差异最小化。基于局部统计特征的方法，根据图像的局部变化与对比，计算单个像素或者局部窗口的阈值。通常对于不同的局部窗口，会得到不同的阈值，属于一种局部阈值化算法。除了上述分类清楚的算法之外，阈值化处理还包括一些其它类型的算法，包括基于模糊集的算法，基于进化的算法，混合优化算法。但是，这几种算法的基本理念与前述几类算法在本质上是一致的，当图像中存在严重噪声污染时，这些算法通常无法得到满意的阈值化结果。

此外，在进行基于血管骨架和边缘提取的血管检测算法，以及血管轮廓提取算法等之前，都需要对血管图像进行阈值化处理，从而提取出血管图像部分。在应用某些血管分割方法之前，也需要进行一个或多个阈值化处理，例如Socher’s方法，一种使用多层分类器的边缘空间学习方法，必须在第一层找到边缘，并且边缘检测的精度取决于应用阈值化算法对污染区域进行消除噪声的结果的质量。因此，阈值化处理是血管检测的一项重要手段。如何提出新的阈值化处理方法，能够有效地去除血管中的噪声，是医学图像领域中的一项重要需求。

血管图像中存在噪声污染问题，全装显微成像和双光子成像能够达到微米或者亚微米的分辨率，能够完成小血管成像，例如毛细血管的成像。毛细血管的直径平均而言在8微米左右，但有可能低至3微米。然而，毛细血管丛通常是缠绕在一起的不规则血管网络，要从噪声背景下自动提取其中的血管是一个 非常具有挑战性的问题。

血管造影图像中，有个特别的问题，即污染区域。血管造影成像通常需要向血管中注入荧光染料，如果血管壁存在裂口，荧光染料会流入周围的组织中。即使在不存在裂口的情况下，长时间的成像，特别是活体成像，都会造成荧光染料向周围组织的渗透。这些渗透的荧光染料，会在血管图像中造成污染区域，类似云团状噪声。这些区域的尺寸通常较血管的尺度大很多。例如，图1a中，血管图的中上位置处的云团状污染噪声区域可以达到130×50μm²。此外，云团状噪声中的像素的亮度在70到165之间(红色通道)，血管像素的亮度则在80到200之间(红色通道)。可以看出，某些噪声像素的亮度值要高于某些血管像素的亮度值。这个问题给准确的血管检测提取带来了困难。然而，现有的阈值化血管提取方法不能很好得解决该问题。其结果要么是真实的血管中的部分像素被错误得剔除了，要么是噪声中的部分像素被错误得保留了下来。不准确的血管提取结果，将给后续的血管图像定量分析(例如，血管长度计算、血管面积或者体积计算等)引入误差。

为了更好得说明血管图像阈值化中存在的问题，以图1a中的血管图像为例，对其用典型阈值化算法处理的结果进行说明。图1a中可以看到在中上部存在明显的云状噪声污染区域。首先应用Kittler’s方法对转化为灰度图后的图1a进行阈值化处理，提取其中的血管，其结果在图1b中给出。之所以选择Kittler’s方法，在于该方法在不同情况下的平均性能要优于其它多种阈值化方法。从图1b中的结果可以看出，虽然大多数血管像素都保留了下来，但大量的云状噪声也都没有去除。此外，为了说明错误得去除亮度较低的小血管的问题，采用最大熵方法对图1a中血管图像进行阈值化处理，其结果在图1c中给出。从其结果中可以看出，云状噪声都消除了，但同时大量的血管也都被错误的去除了。

上述各种问题并不只存在于毛细血管图像中。例如，相似的问题还出现在 虹膜图像、纱槌图像、航拍的道路图像，以及许多其他存在非均匀照明或者引起噪声的物理原因的图像中。因此，有必要开发一种自适应的阈值化处理方法，能够区分有实际物理意义的模式和噪声模式，从而合理的处理污染区域等问题。

发明内容

本发明的目的在于提供一种基于随机方向直方图比的双阈值的血管图像处理方法，该方法基于随机探针检测血管图像子窗口内的图像性质指标，鉴别血管图像中的噪声污染区与血管区域，通过该方法获能够有效得去除血管图像中的噪声。

本发明是通过以下技术方案来实现：

一种基于随机方向直方图比的双阈值的血管图像处理方法，包括以下步骤：

步骤一，输入采集的血管图像，进行背景清理；

步骤二，针对整幅血管图像，应用多级Otsu算法获取两个阈值，分别记作t1和t2，且t1>t2；

步骤三，将步骤一背景清理后的血管图像划分为w_max个子图像，在每个子图像窗口内随机放置n个探针，对n个探针在各个方向上的探测结果进行累加，获得该子图像的方向直方图，然后计算该方向直方图中最大值与最小值的比，得到该子图像的方向直方图比DHR；以此计算得到整幅血管图像中各个子图像的方向直方图比；

步骤四，若某个子图像的方向直方图比DHR大于2，则该子图像为背景干净的血管图像，应用步骤二中所得的低阈值t2对图像进行阈值化处理；

反之则为噪声污染区域，应用步骤二中的高阈值t1对图像进行阈值化处理；

对所有子图像的处理结果汇总，得到了最终的血管图像处理结果。

步骤一对采集的血管图像进行背景清理，具体操作为：对采集的血管图像统计其中最亮像素的亮度值，将该血管图像中所有亮度值低于最亮像素亮度值的5％的像素置为零，其他像素则保留原有的像素值。

步骤二所述的多级Otsu算法采用三级Otsu算法。

步骤三所述的探针，随机设置在子图像窗口内非零像素的位置。

每个探针检测其周围8个或者16个方向上的非零像素的数目，并计算出8维或者16维的方向直方图，在每个方向上，如果探测到的像素值大于零，则在方向直方图的对应维度上加1，反之则不加；对n个探针在各个方向上的探测结果进行累加，共同组成一个方向直方图，该直方图为相应子图像的方向直方图。

与现有技术相比，本发明具有以下有益的技术效果：

本发明的基于随机方向直方图比的双阈值的血管图像处理方法，根据血管的几何分布特点，提出了一种基于随机探针检测血管图像子窗口的图像性质的指标，即方向直方图比，根据方向直方图比的值大小，能够鉴别血管图像中的噪声污染区域与血管区域，对噪声污染区域采用高阈值t1进行阈值化处理，从而能够尽量的消除噪声像素，而对血管区域则采用低阈值t2进行阈值化处理，从而尽量的保留血管像素。通过本发明方法能够有效得去除血管图像中的噪声，且无需复杂的数学模型，计算速度快，去噪能力强，血管图像处理效果好。

附图说明

图1a为双光子血管图像中的云状噪声污染问题结果图；

图1b为Kittler’s方法阈值化处理结果图；

图1c为最大熵方法阈值化处理结果图；

图2为8方向探针统计方向直方图说明(图中的阴影部分为管状物，像素

i是随机选取的管状物上的某个像素，像素j和k是像素i在4和6方向的d-邻居像素)；

图3a为d-邻居像素二维示意图(像素i处的随机探针，j像素为其3-邻居像素)；

图3b为d-邻居像素三维示意图(中心浅色像素表示探针位置，周围的黑色像素表示其1-邻居像素，为显示清晰起见，只标出了六个1-邻居像素)；

图4为基于随机方向直方图比的双阈值血管图像处理方法流程图；

图5为方向直方图生成步骤流程图；

图6为高质量双光子血管图像处理结果比较图；

其中，(a)为原始血管图像图；(b)为本发明算法血管图像处理结果图；(c)为Otsu算法血管图像处理结果图；(d)为Huang算法血管图像处理结果图；(e)为Percentile算法血管图像处理结果图；(f)为Triangle算法血管图像处理结果图；(g)为最大熵算法血管图像处理结果图；(h)为Kittler算法血管图像处理结果图；

图7为中度污染的双光子图像血管像处理结果比较图；

其中，(a)为原始血管图像图；(b)为本发明算法血管图像处理结果图；(c)为Otsu算法血管图像处理结果图；(d)为Huang算法血管图像处理结果图；(e)为Percentile算法血管图像处理结果图；(f)为Triangle算法血管图像处理结果图；(g)为最大熵算法血管图像处理结果图；(h)为Kittler算法血管图像处理结果图；

图8为重度污染双光子图像血管提取结果比较图；

其中，(a)为原始血管图像图；(b)为本发明算法血管图像处理结果图；(c)为Otsu算法血管图像处理结果图；(d)为Huang算法血管图像处理结果图；(e)为Percentile算法血管图像处理结果图；(f)为Triangle算法血管图像处理结果图；(g)为最大熵算法血管图像处理结果图；(h)为Kittler算法血管图像处理结果图；

图9a为三维血管图像原图；

图9b为本发明方法的三维血管图像血管图像处理结果图。

具体实施方式

下面结合具体的实施例对本发明做进一步的详细说明，所述是对本发明的解释而不是限定。

为了解决血管图像中去除噪声污染的问题，本发明基于血管几何分布的特点，即血管像素在某个方向上的分布占优，提出了一种基于随机探针检测血管图像子窗口的图像性质的指标，即方向直方图比，根据方向直方图比的值大小，能够鉴别血管图像中的噪声污染区域与血管区域，对噪声污染区域采用高阈值进行阈值化处理，从而能够尽量的消除噪声像素，而对血管区域则采用低阈值进行阈值化处理，从而尽量的保留血管像素，其中，高、低两个阈值则通过三级Otsu算法获得，通过该方法获得的血管检测结果，能够有效得去除血管图像中的噪声。

本发明的方法称为基于随机方向直方图比的双阈值血管检测方法，简记为DHR-RP(Directional Histogram Ratio at Random Probes)方法。DHR-RP方法基于蒙特卡洛思想，采用先随机抽样后统计分析的方式，提取血管图像中的局部区域特征。进一步，通过区域特征的差异(血管存在区域和噪声污染区域)来选择不同的阈值。

本发明采用以下技术方案予以实现，一种基于随机方向直方图比的双阈值的血管图像处理方法，其特征在于，包括以下步骤：

步骤一，对图像背景进行清理，统计图像中最亮的像素的亮度值，将血管图像中所有亮度值低于最亮像素亮度值的5％的像素置为零，其他像素则保留原有的像素值，由于探针的位置会随机选择像素值非零(亮像素)的位置，该处理使得随机选取的探针能够位于真正的亮像素位置上，避免由于背景位置成像的像素值为接近零的非零值，而导致放置大量无用探针，降低计算效率；

步骤二，针对整幅血管图像，应用多级Otsu算法计算高、低两个阈值，以用于后续的双阈值化处理，本发明选取三级Otsu算法获取两个阈值t1和t2，分别设为高阈值t1和低两阈值t2；

步骤三，将血管图像划分为w_max个子图像，在每个子图像窗口内随机放置n个探针，其位置在子图像窗口内的非零像素位置随机选取，在每个探针的位置检测其周围8个或者16个方向上检测非零像素的数目，从而计算8维或者16维的方向直方图，在每个方向上，如果探测到的像素值大于零，则在方向直方图的对应维度上加一，反之则不加，n个探针在各个方向上的探测结果进行累加，即对相应方向统计求和，共同组成一个方向直方图，该直方图即相应子图像的方向直方图，然后计算方向直方图中最大值与最小值的比，该比值定义为方向直方图比(Directional Histogram Ratio,简记为DHR)，按照该方法计算各个子图像的方向直方图比；

步骤四，对子图像进行双阈值化处理，子图像的方向直方图比表示了相应子图像的性质，经理论证明，如果方向直方图比大于2，则子图像为背景干净的血管图像，反之则为噪声污染区域，因此，双阈值化处理时，如果子图像中方向直方图比DHR>2，应用步骤二中所得的低阈值t2对图像进行阈值化处理，反之，则应用步骤二中的高阈值t1对图像进行阈值化处理，所有子图像的处理结果组成了最终的血管检测结果。

基于随机方向直方图比的双阈值血管图像处理方法(DHR-RP)包含多个步骤，主要步骤包括背景清理、子窗口分割、计算直方图比、双值滤波等。

首先，详细介绍DHR-RP方法的基本概念、思路。

相对比而言，噪声像素通常随机得分布于整个图像中，而某个物体对应的像素通常是集中分布在有一定几何边界的范围内，存在某种与其形状对应的特 定模式，例如方形、圆形、关节连接等。DHR-RP方法正是基于血管特殊形状(管状)的几何分布而建立的。

具体来讲，作为一种管状物，血管在轴向的长度远远大于径向的尺寸。以此为先验知识，如果放置一个星形的探针在某个血管像素上，如图2所示，在每个方向上对亮的像素(背景通常为暗的像素)进行采样统计。一个合理的猜测是，在血管的轴向统计所得的亮像素的个数要远远大于在轴向统计所得的个数。需要注意的是，星形探针采样的臂长需要大于血管的直径，从而能够获取到有效的形状分布信息。然而，当星形探针在云状噪声中采样时，由于噪声分布的随机性，其在各个方向上统计所得的亮像素的个数会趋于相近。因此，星形探针统计所得到各个方向上的亮像素的个数的分布，可以作为一个判断某个像素是血管像素，还是随机噪声的指标。在这种假设的前提下，提出了DHR-RP算法。

下面给出d-邻居像素的定义：

定义1(d-邻居像素)：设i和j是一帧图像中的两个像素，且i≠j。如果在沿i→j的直线方向上，i和j之间的距离为d，则称像素j是像素i在i→j方向的d-邻居像素。

像素i所有方向上的d-邻居像素组成的像素集合，表示为  $N_i^d=\{N_{\mathrm{ij}}^d,j=\mathrm{1,2},...,h\}(j=\mathrm{1,2},...,h).$ 其中表示像素i在i→j方向上的d-邻居像素，h是所统计的方向的数量。

在定义方向直方图之前，需要先给出随机探针的概念。简单来讲，每个随机探针都随机得放置在一个具有大于零的亮度值的像素处，并利用星形探针来统计各个方向的亮度大于零的像素的数量。因为一般而言，血管像素或者噪声通常具有大于零的亮度值，而背景通常是暗色的，在经过简单的预处理之后，亮度值变为零。接下来给出方向直方图的定义。

定义2(方向直方图)对第i(i＝1,2,...,n)个像素，以表示其第j(j＝1,2,...,h)个方向的d-邻居像素。设具有h维的直方图，其第j维统计所有随机探针在第j方向上具有大于零的亮度的邻居像素的数量，则该直方图称为方向直方图，表示为其中

$H_j^d=\underset{i=1}{\overset{n}{\mathrm{\Σ}}}1(I_{\mathrm{ij}}^d>0),$    公式(一)

其中，1(·)是指示函数。

为了更好的说明方向直方图的概念，图3a和图3b给出了d-邻居像素和随机探针在二维和三维两种情况下的示意图。在图3a中，中心的浅色像素i是随机选择的像素位置，其最外圈的8个像素分别是像素i的d-邻居像素，d＝3。在图3b中，中心的浅色圆为随机探针所在位置，其周围的深色圆为体素(立体像素)i的d-邻居像素，d＝1。此外，图3a和图3b中的中心像素的随机探针，作为二维和三维的情况，可以分别选择8个或者26个直线方向，此时对应的方向向量。如果第j个方向的d-邻居像素的亮度大于零，则在方向向量对应的第j维加1。在二维图像中j＝1,2,…,8，在三维图象中j＝1,2,…,26。

其次，详细解释方向直方图比的概念。

为了描述各个方向上d-邻居像素的数量差，定义了方向直方图比。方向直方图比的值将作为选择阈值的标准。

定义3(方向直方图比)方向直方图向量中的最大值和最小值之比，称为方向直方图比，表示为

再次，详细推导应用方向直方图比进行血管检测阈值选择的方法。

为了证明方向直方图比可以作为阈值选择的标准，需要首先讨论管状物的 几何特征。简单起见，我们首先考虑图2中的二维图像的例子。如图2所示，设像素i是某个管状物O(阴影部分)中的一个像素，表示为i∈O，像素j和k分别为像素i在i→j方向和i→k方向的d-邻居像素。在图2中，i→j方向选择为血管的轴向，而i→k方向为非轴向。如果d的值大于管状物的直径，则在除了轴向方向的其他方向上(例如，i→k方向)，d-邻居像素是处于管状物的外部，包括图2中的1、2、3、5、6、和7方向。因此，对于图2中的管状物O，对任意像素i∈O，总是存在不同方向的d-邻居像素j和k，使得j∈O且 对于形状呈管状的血管而言，该结论至少在血管的局部区域中成立。进一步而言，由于像素i在管状物之内，因此具有大于零的亮度值，在方向直方图的相应维数上应当加1；当像素k落在背景区域时，该像素的亮度为零，则其对方向直方图没有贡献。因此可以推断，一个存在管状物的区域所统计的方向直方图中(d-邻居像素)，对应管状物轴向的相应的维数的值将远远大于非轴向所对应的各维数值。因此，其对应的方向直方图比在理想情况下(如图2)趋于无穷大，因为方向直方图中最小值为零，即min(H^d)＝0。然而，在真实的血管图像中，存在噪声，而且通常有多条血管，因此方向直方图的最小值min(H^d)不精确为零，而是某个很小的值。但此时，方向直方图比的值仍然比较大。而对于污染区域，噪声像素(亮度非零)的分布相对均匀，因此计算所得的方向直方图的值在各个维度上的值将十分接近，从而会得到一个较小的方向直方图比。

综上所述，方向直方图比DHR可以作为管状物是否存在的判断指标。为了从理论上论证该结论，这里考虑血管图像中的两种情景。情景1：图像中存在多条管状物，背景干净；情景2：存在很少或者没有管状物，背景污染。以下分别就这两种情景做出证明。

i是随机选取的管状物上的某个像素，像素j和k是像素i在4和6方向的d-邻居像素)

命题1设图像中的物体都是管状物，且其空间分布不与统计方向直方图的所有方向具有各向同性，以E(Dr)表示方向直方图比的期望，则对于情景1存在E(Dr)＞＞2，对于情景2存在1≤E(Dr)<2。

证明：设m表示图像中的总像素数，其中m₁个像素来自图像中的物体(像素值大于零)，m₂个像素来自噪声(像素值大于零)，m₃个像素来自背景(像素亮度为零)。分别以M₁和M₂表示m₁个物体像素和m₂个噪声像素对应的像素集合。根据随机探针的选择方法，随机探针的位置将只会在集合M₁∪M₂中。第i(i＝1,2,...,n)个随机探针来自于某个物体的概率为

$p(i\&Element;M_1)=\frac{m_1}{m_1+m_2},$    公式(二)

而来自于噪声的概率为：

$p(i\&Element;M_2)=\frac{m_2}{m_1+m_2}.$    公式(三)

设表示随机探针i的d-邻居像素。为简化分析，不考虑随机探针放置像素i和之间的相关性。此时，像素来自某个物体的概率满足如下不等式

$0\&le;p(N_{\mathrm{ij}}^d\&Element;M_1|i\&Element;M_1\&cup;M_2)\&le;\frac{m_1}{m},j=\mathrm{1,2},...,h.$    公式(四)

其中，不等式中的下界对应图2中d-邻居像素位于i→k方向的情况，上界对应i→j方向的情况。由于噪声像素在整个图像中随机分布，因此像素属于噪声的概率为

$p(N_{\mathrm{ij}}^d\&Element;M_2|i{\&Element;M}_1\&cup;M_2)=\frac{m_2}{m}.$    公式(五)

如果d-邻居像素的亮度值大于零，则在方向直方图中所对应的第j维加1。设表示一个事件，即亮度为正且在方向直方图中的第j为得以计 数。根据公式(四)和(五)可以得出事件发生的概率为

$\begin{array}{c}p\left(A_{\mathrm{ij}}^d\right)=p(N_{\mathrm{ij}}^d\&Element;M_1|i\&Element;M_1\&cup;M_2)\\ +p(N_{\mathrm{ij}}^d\&Element;M_2|i\&Element;M_1\&cup;M_2)\\ \&DoubleRightArrow;\frac{m_2}{m}\&le;p\left(A_{\mathrm{ij}}^d\right)\&le;\frac{m_1+m_2}{m}\end{array}.$    公式(六)

在方向直方图的第j维做出贡献的d-邻居像素总数的期望满足

$\begin{array}{c}E\left(H_j^d\right)=\underset{i=1}{\overset{n}{\mathrm{\&Sigma;}}}p\left(A_{\mathrm{ij}}^d\right)\&CenterDot;1\\ \&DoubleRightArrow;n\&CenterDot;\frac{m_2}{m}\&le;E\left(H_j^d\right)\&le;n\&CenterDot;\frac{m_1+m_2}{m}\end{array}.$    公式(七)

根据公式公式(七)以及方向直方图比的定义可知，方向直方图比的期望如公式(八)所示：

$E\left(\mathrm{Dr}\right)=\frac{\underset{j}{\max }\left(E\left(H_j^d\right)\right)}{\underset{j}{\min }\left(E\left(H_j^d\right)\right)}=1+\frac{m_1}{m_2}.$    公式(八)

在情景1中，存在m₁＞＞m₂，而且对应的m₂很小，接近为零。此时，由公式(八)可得

$E\left(\mathrm{Dr}\right)=1+\frac{m_1}{m_2}>>2.$    公式(九)

在情景2中，存在m₂＞＞m₁，因此有

$1\&le;(E\left(\mathrm{Dr}\right)=1+\frac{m_1}{m_2})<2.$    公式(十)

因此，在这两种情况下(血管区域和污染区域)，方向直方图比位于不同的取值区间，可以应用方向直方图比判断血管存在与否。

最后，介绍本发明的具体实施示例，基于随机方向直方图比的双阈值图像处理方法的整体步骤如图4所示，计算方向直方图的步骤如图5所示。

步骤一，图像背景清理，统计图像中最亮的像素的亮度值，将血管图像中所有亮度值低于最亮像素亮度值的5％的像素置为零，其他像素则保留原有的像素值；

步骤二，针对整幅血管图像，应用多级Otsu算法计算高、低两个阈值，以用于后续的双阈值化处理，本发明选取三级Otsu算法获取两个阈值分别记作t1和t2，且t1>t2；分别记为高阈值t1、低阈值t2；

步骤三，由于DHR-RP算法是一种局部方法，其执行过程中需要将图像划分为多个子窗口中，在每个子窗口中，散布随机探针，统计相应的方向直方图。在具体实施中，二维图像和三维图像中的子窗口分别采用长方形和长方体。在每个子窗口内，方向直方图通过图5中的各个步骤产生，其中rp随机探针数量的计数器，n是预先指定的某个子窗口图像中的随机探针数量的最大值，n的值默认为100，在子窗口中随机放置一个星形探针之后，沿如图2中所示的各个方向搜索其d-邻居像素，如果d-邻居像素的亮度大于零，则在方向直方图中的相应维加1，反之则维持方向直方图不变，最后，计算每个子图像的方向直方图比；

步骤四，对子图像进行双阈值化处理，如果子图像中的方向直方图比DHR大于2，应用步骤二中所得的低阈值t2对图像进行阈值化处理，反之，则应用步骤二中的高阈值t1对图像进行阈值化处理，所有子图像的处理结果组成了最终的血管检测结果，此外，在三维图像中，双阈值需要针每幅二维子帧图像计算。

图6、图7和图8分别给出了应用本发明方法时的，高质量双光子血管图像、中度污染双光子血管图像和重度污染双光子血管图像的血管图像处理结果，以及与其它方法的比较。三组图中，(a)为原始血管图像图；(b)为本发明算法血管图像处理结果图；(c)为Otsu算法血管图像处理结果图；(d)为 Huang算法血管图像处理结果图；(e)为Percentile算法血管图像处理结果图；(f)为Triangle算法血管图像处理结果图；(g)为最大熵算法血管图像处理结果图；(h)为Kittler算法血管图像处理结果图；图9给出了一组应用本发明方法的三维血管图像的检测结果。这些结果显示了本发明方法在血管检测中优秀的去除噪声能力，以及能够获得清晰高质量的图像处理结果。

Claims (5)

1.一种基于随机方向直方图比的双阈值的血管图像处理方法，其特征在于，包括以下步骤：

步骤一，输入采集的血管图像，进行背景清理；

步骤二，针对整幅血管图像，应用多级Otsu算法获取两个阈值，分别记作t1和t2，且t1＞t2；

步骤三，将步骤一背景清理后的血管图像划分为w_max个子图像，在每个子图像窗口内随机放置n个探针，对n个探针在各个方向上的探测结果进行累加，获得该子图像的方向直方图，然后计算该方向直方图中最大值与最小值的比，得到该子图像的方向直方图比DHR；以此计算得到整幅血管图像中各个子图像的方向直方图比；

步骤四，若某个子图像的方向直方图比DHR大于2，则该子图像为背景干净的血管图像，应用步骤二中所得的低阈值t2对图像进行阈值化处理；

反之则为噪声污染区域，应用步骤二中的高阈值t1对图像进行阈值化处理；

对所有子图像的处理结果汇总，得到了最终的血管图像处理结果。

2.根据权利要求1所述的一种基于随机方向直方图比的双阈值的血管图像处理方法，其特征在于，步骤一对采集的血管图像进行背景清理，具体操作为：对采集的血管图像统计其中最亮像素的亮度值，将该血管图像中所有亮度值低于最亮像素亮度值的5％的像素置为零，其他像素则保留原有的像素值。

3.根据权利要求1所述的一种基于随机方向直方图比的双阈值的血管图像处理方法，其特征在于，步骤二所述的多级Otsu算法采用三级Otsu算法。

4.根据权利要求1所述的一种基于随机方向直方图比的双阈值的血管图像处理方法，其特征在于，步骤三所述的探针，随机设置在子图像窗口内非零像素的位置。

5.根据权利要求4所述的一种基于随机方向直方图比的双阈值的血管图像处理方法，其特征在于，每个探针检测其周围8个或者16个方向上的非零像素的数目，并计算出8维或者16维的方向直方图，在每个方向上，如果探测到的像素值大于零，则在方向直方图的对应维度上加1，反之则不加；对n个探针在各个方向上的探测结果进行累加，共同组成一个方向直方图，该直方图为相应子图像的方向直方图。