CN101882222B - 基于基元结构定义与区域生长技术的虹膜分区及阳光放射沟提取方法 - Google Patents

Abstract

本发明提出了一种基于基元结构定义与区域生长技术的虹膜分区及阳光放射沟提取方法，包括一种融合中西医虹膜学虹膜分区的新虹膜分区图谱(CADIC)；并发明了应用此图谱所必需的基于基元定义的区域生长技术用于典型虹膜特征——阳光放射沟的提取，并发明了应用此图谱所必需的基于虹膜图像动态分区技术和虹膜图谱的自适应覆盖方法。本发明在保证了肌体各部脏器在虹膜上的准确定位同时，突出了虹膜图谱的实用性，实现了医学意义下完整虹膜特征的动态分割技术及其在虹膜图谱上的定位技术。通过实验证明，此技术能够大大提高了虹膜识别和诊病可靠性。

Description

基于基元结构定义与区域生长技术的虹膜分区及阳光放射沟提取方法

(一)技术领域

本发明涉及计算机图像处理技术领域，具体涉及一种利用虹膜进行疾病诊断的技术。

(二)背景技术

目前，随着社会的进步和人民生活水平的日益提高，人们越来越关注自身的健康，迫切需要一种方便可靠、无创伤非接触式的自我诊断手段。虹膜学与虹膜诊断正是这样一种现代新型的诊病方法，它具有非接触性、无伤痛性、易于普及的特点，这也是未来医学与预防医学诊断手段的重要发展趋势。

中西医虹膜学认为：人眼虹膜是人体中与全身各脏器密切相关的器官之一，各部脏器在虹膜上有确切的映射区，其变化情况直接反映人体内部各脏器的健康状况。虹膜异常主要体现在各个脏器映射区内，色调的变化、局部色斑的出现、虹膜纤维破洞的出现或虹膜纤维疏密程度的变化等方面。通过观察虹膜不同区位的异常变化，可以对某些器官的疾病做出诊断，如：虹膜上方脑映射区内某种亮点表示脑神经有异常、肠胃映射区瞳孔周围出现凹点说明有溃疡病等。

利用虹膜诊病的关键所在是定位各部脏器在虹膜上的映射区。为此中西医虹膜学均给出了虹膜图谱。从西方虹膜学之父匈牙利医生Ignatz Von Peczley到美国虹膜大师BernardJensen，西方虹膜学从实验观察和统计分析出发，其图谱定位准确，所以提出的诊断标准比较明确和稳定。但西医对虹膜分区过于复杂，时常将一个完整病症被分割在几个区域中，没有摆脱“人是无数零件组成的机器”的观点，因此在临床上无法把这些活体器官之间的联系统一起来；另一方面西医分区在利用计算机辅助诊断方面的可操作性较差，因为目前的虹膜图像采集设备和虹膜图像处理算法上还不能达到直接利用西医图谱诊病所要求的精度。

我国中医“眼诊学”认为：眼虽是一个局部器官，但它与全身的脏腑、经络共同构成了一个有机的整体：在生理上相互协调，在病理上相互影响，有着非常密切的关系。所以中医虹膜图谱在整体上给出了映射区在虹膜上的大致位置。因此造成了中医图谱精度较差，并且目前还没有一个统一的、严格的和量化的标准。

此外，中西医图谱还存在一个共同的问题：病灶出现的位置会跨越相邻的两个分区，出现这种情况包括两方面的情况：一、相关脏器(如生殖泌尿系统中的各脏器)出现在相邻的映射区中，从疾病的发生、发展来考察，他们是相互关联、相互影响的，其病灶也有可能扩大到其它映射区中；另一方面，某个脏器病变较重(如重症冠心病，其心脏映射区病理变化较大)。

(三)发明内容

本发明的目的在于提供一种有机结合中、西医的虹膜学成果---依据西医虹膜学中较准确定位肌体脏器与虹膜映射区域的关系、根据中医理论对虹膜映射区形成的阐述，给出基于基元结构定义与区域生长技术的虹膜分区及阳光放射沟提取方法。

本发明的目的是这样实现的：

1、基于中西医虹膜学的计算机辅助诊断虹膜图谱(CADIC：computer aided diagnosis irischar)

本发明的虹膜图谱定义如下：对于一幅虹膜图像，经过前期预处理与虹膜内外边缘精确定位后，将虹膜图像置于极坐标系(θ，r)中，把虹膜分为以瞳孔圆心为极点、以瞳孔边缘为起始边界的五个环、12个区的虹膜分区图谱。

设一幅虹膜图像用A表示，R_T为瞳孔半径，R_I为虹膜半径，R_D为图像中卷缩轮的理论半径，R′_D为卷缩轮的实际半径。极坐标的极轴方向取为时钟的零点方向，极角按照顺时针方向定义。

将A划分为五个环12个区后，其中的任意一个区用A_i ^k(θ，r)来表示。其中：k＝1，2，...，12；i＝1，2，...，5，每个分区的定义如下：

A_i ^k(θ，r)： $\frac{(k-1)\mathrm{\&pi;}}{6}<\mathrm{\&theta;}\&le;\frac{\mathrm{k\&pi;}}{6}$

$\left\{\begin{array}{c}\mathrm{if}:i=1;\mathrm{then}:r=R_T\\ \mathrm{if}:i=2;\mathrm{then}:R_T<r\&le;R_T+\frac{2}{5}(R_I-R_T)\\ \mathrm{if}:i=\mathrm{3,4,5};\mathrm{then}:R_T+(i-1)\frac{(R_I-R_T)}{5}<r\&le;R_T+i\frac{(R_I-R_T)}{5}\end{array}\right.$

虹膜分区是按时钟顺时针方向从12点标记开始每30度为一个区；在虹膜环形的半径方向上5等分，构成5个环。根据虹膜学理论，虹膜中卷缩轮大小占虹膜的三分之一左右，所以最里边的两个环被合并表达卷缩轮部分。再根据虹膜的另外三个环和分区，将虹膜分成各个诊断区域。本发明根据制定的分区规则和医学虹膜图谱进行了比照，得到了一一对应关系，表1-8给出了左右眼虹膜CADIC图谱中各个分区对应的肌体脏器在虹膜上的映射部位。

2、基于基元定义和区域生长技术的虹膜处理方法

本发明提出的动态自适应虹膜分区(图谱覆盖)方法针对虹膜图像纹理复杂、区域变化不明显特性，在基于串行区域生长与合并技术基础上，通过对特征先验知识，定义区域生长与合并准则，以“基元”为区域生长点，具体地方案包含如下步骤：

1.虹膜图像的采集；

2.虹膜图像的预处理；

3.虹膜定位；

4.CADIC虹膜图谱覆盖；

5.基于基元定义和区域生长的阳光放射沟提取；

6.虹膜自适应分区。

各个处理步骤的细节如下：

1.虹膜图像的采集

采样方法是对视频图像抓拍；采样者头部固定并正视；采样时所采用的光源产生的光斑打在瞳孔和虹膜的交界处；采样图像为768*576的24位位图。每个采样者采3到5组虹膜图像。

2.虹膜图像的预处理

在对虹膜图像进行分区分析之前，要进行预处理工作。虹膜图像的预处理主要包括彩色图像灰度转换、噪声处理、几何校正、图像增强。

2.1)图像灰度化

对采集到的彩色图像，应用如下公式将其变换为灰度图像。

$\left[\begin{array}{c}Y\\ I\\ Q\end{array}\right]=\left[\begin{array}{ccc}0.299& 0.587& 0.114\\ 0.596& -0.274& -0.322\\ 0.211& -0.523& 0.312\end{array}\right]\left[\begin{array}{c}R\\ G\\ B\end{array}\right]---(2-1)$

Y＝0.299×R+0.587×G+0.114×B    (2-2)

式中R——真彩图中每个像素的红分量；

G——真彩图中每个像素的绿分量；

B——真彩图中每个像素的蓝分量；

Y——该像素所对应的灰度

2.2)虹膜图像灰度归一化

通过将原始虹膜图像中感兴趣的区域映射到相同的尺寸和相同的对应区域中进行虹膜归一化，从而消除平移、旋转对纹理检测的影响。

2.3)光斑填充

采用数字图像的自适应插值法来对光斑进行填充。零阶插值算法，即令输出像素的灰度值等于离它所映射的位置最近的输入像素的灰度值。

3.虹膜定位

除去那些于虹膜无关的信息(如眼睑Eyelid，睫毛Eyelash，瞳孔Pupil等)，也就是从虹膜图像中定位出虹膜的内边界和外边界。

3.1)虹膜粗定位

通过对虹膜图像的分析，得出以下的先验知识：瞳孔的直径最小也要大于120个像素；瞳孔是眼睛图像中面积最大的低灰度区域，而其圆心位置应该在上下左右灰度和值最小的区域附近，利用这些先验知识，可以减小搜索的空间提高精度。本文采用传统的灰度投影法对虹膜内、外边缘粗定位，得到虹膜内径r、外径r1和瞳孔中心(x_c，y_c)；

3.2)虹膜精定位

本发明采用Daugman算法来对虹膜进行精定位，其核心在于利用一个有效的积分微分算子来计算圆参数。虹膜具有良好的环状特性，所以可采用圆形探测算子对虹膜边界进行精定位。

在虹膜粗定位瞳孔中心(x_c，y_c)和半径r的基础上，采用圆探测算子实现虹膜精定位；所述虹膜精定位方法为基于圆探测算子的定位方法：

a、对于中心矩阵内的每一个元素，分别计算半径为r同一圆周上的像素点灰度和sum，其中中心矩阵的范围是(x_c-5，y_c-5)，半径变化范围为r-5～r+5；

b、对于r的变化值d，计算出对应的高斯函数值：

$g\left(d\right)=(1/(0.1\&times;\left(\sqrt{2\mathrm{\&pi;}}\right)))e^{((-{(d-6)}^2)/(2\&times;{0.1}^2))};$

c、计算灰度梯度值向量：det＝|g(d+1)-g(d)|×sum；

d、依次计算所有元素对应的det向量值；

e、在所有的det向量中取元素最大值，即为精定位的结果。

4.CADIC虹膜图谱覆盖

在虹膜内外边缘精确定位后，将虹膜部分所在区域的直角坐标(x，y)形式转换为极坐标(r，θ)形式。根据虹膜图谱和r、θ的值来将虹膜部分分成若干个区域。这是第一次分区，也就是理论分区。

5.基于基元定义和区域生长的阳光放射沟提取

本发明选取虹膜图像中疾病形成的阳光放射沟为处理对象，首先定义纹理基元并讨论基于基元的纹理表达方法。

设：纹理基元为h(x，y)，排列规则为r(x，y)，则纹理t(x，y)定义为

$t(x,y)=h(x,y)\&CircleTimes;r(x,y)---(4-1)$

设x_m，y_m为脉冲函数的位置坐标，则

r(x，y)＝∑δ(x-x_m，y-y_m)    (4-2)

根据卷积定理，在频域有：

T(u，v)＝H(u，v)R(u，v)      (4-3)

R(u，v)＝T(u，v)H(u，v)^-1    (4-4)

那么对于给定纹理基元h(x，y)的描述，即可得到反卷积滤波器H(u，v)^-1。将这个滤波器用于待处理纹理图像，每个脉冲都在纹理基元的中心。对整幅图像中不同纹理基元的分布进行统计可获得图像的全面纹理信息。

通过纹理基元确定对了对待分割的区域找到一个种子子图---即基元脉冲中心，然后将种子像素周围邻域中与种子像素有相同或相似性质的像素(通过定义的相同或相似规则)合并到种子像素所在区域中。在将新像素作为种子继续扩展，直到再没有满足条件的像素可被包括，这样的一个区域就长成。

对于一个给定的目标图像X和一个结构基元S，以符合S的种子区域为中心在图像上进行迭代生长。

对于种子

的每一个邻近像素x，在生长过程中应始终保持：

if：S(x)＝TRUE            (4-5)

$\mathrm{then}:X\&CirclePlus;S\&LeftArrow;X\&CirclePlus;S\&cup;\left\{x\right\}---(4-6)$

即，如果区域邻近像素x能符合结构基元定义，则将像素x合并到特征区域，使之成为新种子的一部分。

如果对于种子

的所有邻近像素x，均有：

S(x)＝FALSE               (4-7)

则生长停止。

对于阳光放射沟特征，通过前期对临床病理的研究分析，从获得的先验统计知识定义像素生长规则如下。

规则1：沟长度为L个像素，5≤L≤R_虹膜；

规则2：沟宽度为W个像素，3≤W≤5；

规则3：沟的延长线与瞳孔交点的切线相垂直；

规则4：根据虹膜图像灰度直方图确定沟区域灰度最优阈值T_G；

规则5：沟的起点为卷缩轮边界点。

根据CADIC虹膜分区图，对每个感兴趣的诊断候选区域，执行如下阳光放射沟检测算法。

1)计算阳光放射沟灰度阈值：

对整个虹膜统计灰度直方图，定义灰度阈值T_G，使得灰度小于T_G的像素数目为整个虹膜像素数目的10％。根据不同的虹膜图像采集设备，此比例数应根据人工标注阳光放射沟情况的先验统计数据做出相应的调整。

2)对任一个候选区Q，执行以下3)～10)。

3)沟起点选取：

对于任意像素x∈Q，如果有符合规则5的像素x的灰度值f(x)≤t_G，则认为x是一个沟起点。把找到的所有沟起点记为d_n，n＝1，2，3...N。

4)如果N＝0，即没有找到任何沟起点，则说明没有阳光放射沟存在，转10)结束本算法。

5)沟起点筛选与合并，形成种子候选：

对沟底点记d_n，n＝1，2，3...N，先后执行生长与合并操作，合并成种子候选d’_n，n＝1，2，3...N’，N’ ≤N。

6)种子候选筛选，形成种子：

对每个种子候选d’_n定义它在平行于虹膜内边缘方向上的像素跨度数为其宽度W(d’_n)，删除所有W(d’_n)＞5的种子候选d’_n。

所有剩余的种子候选d’_n成为种子，记为z_m，m＝1，2，3...M。在不造成混淆的情况下，也用z_m代表本种子的中心像素。

7)如果M＝0，即没有找到任何种子，则说明没有阳光放射沟存在，转10)结束本算法。

8)区域生长，形成阳光放射沟候选：

在每个种子z_m，m＝1，2，3...M周围邻域中搜索同种子像素具有相同或相似性质的像素，合并到种子像素所在区域中；在将新像素作为种子继续迭代扩展，直到再没有满足条件的像素可被包括，这样一个区域就长成了。

8.1)定义每个种子(中心)z_m的邻域为NB(z_m)，NB(z_m)应符合规则1，2和3；

8.2)定义基元结构为S，此处的S定义为典型阳光放射沟的横截线上的像素集合。

8.3)对于任意像素x∈NB(z_m)，根据公式(4-5)、(4-6)所定义的生长规则进行区域生长。

9)阳光放射沟筛选：

对获得的阳光放射沟候选，根据规则1和2进行筛选。

10)输出结果，结束算法。

通过上述方法实验，可以获得满足裂缝基元定义的线条，即为阳光放射沟。

6.虹膜动态自适应分区

具体做法体现在对病理特征(此处为阳光放射沟，也可以推广到其他虹膜特征上)在分区边界上的处理：当分区时病理特征被分割开来，不能被完整地表达在同一个区域中时，进一步从保持病理特征完整性的角度，重新确定该特征所覆盖区域。

6.1)首先，通过估计病理特征的区域进行病理特征筛选：

对病理特征边缘形状的提取采用估算的方法。求出组成特征的像素点在极坐标中的角度和极径的最大和最小值，这样就大概确定出了纹理所覆盖的区域的面积。从前面特征提取的结果可以看到，只有一部分重要的纹理特征(阳光放射沟)具有诊病价值。我们根据所确定的每个特征区域的面积，选取最重要的若干个特征，比如5个；删除其余的不重要特征。

6.2)结合图15，进行CADIC图谱覆盖，方法如前述。

6.3)判断是否有跨区域的病理特征：

对每一个虹膜图谱上的分区注意考察，如果发现其边界同重要特征有交叠，则记录这个交叠。

6.4)进行区域动态调整，完成自适应分区：

如果一个纹理特征覆盖多个区域，进行如下处理：对每个虹膜图谱分区，统计区域外内特征的像素点数目M；之后对其每处交叠边界，考察此区域外该特征的像素点并进行统计数目M’。

如果区域外像素点中只有极少部分M’＜＜M(例如M小于5％的M)，则对这些少数像素点所分布的区域不予考虑；

如果区域外像素点过于重要M’＞M(说明该特征不是主要处于本区域)，则同样不予考虑，并删除本区域中的特征；

此外，说明区域外特征比较重要；则将此边界邻接的区域同本区域合并，实现了重新的区域划分，也就是进行了自适应分区。

本发明提出了一种融合中西医虹膜学虹膜分区的新虹膜分区图谱(CADIC)；并发明了应用此图谱所必需的基于基元定义的区域生长技术用于典型虹膜特征---阳光放射沟的提取，并发明了应用此图谱所必需的基于虹膜图像动态分区技术和虹膜图谱的自适应覆盖方法。本发明在保证了肌体各部脏器在虹膜上的准确定位同时，突出了虹膜图谱的实用性，实现了医学意义下完整虹膜特征的动态分割技术及其在虹膜图谱上的定位技术。通过实验证明，此技术能够大大提高了虹膜识别和诊病可靠性。本发明涉及到的技术领域属于计算机图像处理、模式识别以及中西医结合与现代化领域，有机结合中、西医的虹膜学成果，根据西医虹膜学叫准确定位肌体脏器与病症区域、根据中医理论解释、并处理病症的不同表现，使中医诊疗方法、技术更好地为大众服务，也使中医理论精萃得以传承和发扬，利用计算机强大的信息处理能力，解决中医诊疗在现代化和信息化建设上面临的困难，这是本发明的主要意义。采用本发明，可以在保证了肌体各部脏器在虹膜上的准确定位，同时有效检测虹膜特征，从而大大提高了虹膜识别和诊病的可靠性与实用性。

(四)附图说明

图1为中医虹膜图谱；

图2为Jensen虹膜图谱；

图3-4为CADIC图谱示意图；

图5为具体实施方式流程图；

图6为彩色虹膜图像；

图7为灰度化后的图像；

图8-9为光斑填充前后对比图像；

图10为虹膜粗定位图；

图11为内、外边缘精定位图；

图12为理论分区图；

图13为图像区域分布图；

图14为使用基元定义的裂缝识别结果图；

图15为CADIC图谱覆盖结果图；

图16为虹膜动态自适应分区结果示意图---阳光放射沟的自适应分区结果图。

(五)具体实施方式

下面结合附图对本发明作进一步的说明：

1：基于中西医虹膜学的计算机辅助诊断虹膜图谱(CADIC：computer aided diagnosis irischar)

计算机辅助诊断虹膜图谱分区方法，既要避免西医虹膜图谱的复杂与实用性差的弊端，又要克服中医图谱精度不足和定位不准的弱点；其分区方法应能够保证各大脏器在映射区的完整性和定位的精确性，还要利于计算机实现又方便样本图像的分类存储。

根据图1，图2，本实施例的虹膜图谱定义如下：对于一幅虹膜图像，经过前期预处理与虹膜内外边缘精确定位后，将虹膜图像置于极坐标系(θ，r)中。把虹膜分为以瞳孔圆心为极点、以瞳孔边缘为起始边界的五个环、12个区的虹膜分区图谱，如图3所示。

设一幅虹膜图像用A表示，R_T为瞳孔半径，R_I为虹膜半径，R_D为图像中卷缩轮的理论半径，R′_D为卷缩轮的实际半径。极坐标的极轴方向取为时钟的零点方向，极角按照顺时针方向定义。

将A划分为五个环12个区后，其中的任意一个区用A_i ^k(θ，r)来表示。其中：k＝1，2，...，12；i＝1，2，...，5，每个分区的定义如下：

A_i ^k(θ，r)： $\frac{(k-1)\mathrm{\&pi;}}{6}<\mathrm{\&theta;}\&le;\frac{\mathrm{k\&pi;}}{6}$

$\left\{\begin{array}{c}\mathrm{if}:i=1;\mathrm{then}:r=R_T\\ \mathrm{if}:i=2;\mathrm{then}:R_T<r\&le;R_T+\frac{2}{5}(R_I-R_T)\\ \mathrm{if}:i=\mathrm{3,4,5};\mathrm{then}:R_T+(i-1)\frac{(R_I-R_T)}{5}<r\&le;R_T+i\frac{(R_I-R_T)}{5}\end{array}\right.$

虹膜分区是按时钟顺时针方向从12点标记开始每30度为一个区；在虹膜环形的半径方向上5等分，构成5个环。根据虹膜学理论，虹膜中卷缩轮大小占虹膜的三分之一左右，所以最里边的两个环被合并表达卷缩轮部分。再根据虹膜的另外三个环和分区，将虹膜分成各个诊断区域。我们根据制定的分区规则和医学虹膜图谱进行了比照，得到了一一对应关系，可见，我们所提出的实用化虹膜图谱同传统的中西医虹膜图谱是十分吻合的。见图4，表1-8给出了左右眼虹膜CADIC图谱中各个分区对应的肌体脏器在虹膜上的映射部位。

表1右眼虹膜CADIC图谱

表2右眼虹膜CADIC图谱

表3右眼虹膜CADIC图谱

表4右眼虹膜CADIC图谱

表5左眼虹膜CADIC图谱

表6左眼虹膜CADIC图谱

表7左眼虹膜CADIC图谱

表8左眼虹膜CADIC图谱

2、基于基元定义和区域生长技术的虹膜处理方法

理论分区是通过对正常人虹膜图像学习而制定的，如果人体有病变发生，在对应的虹膜区域就会出现明显的纹理变化。而随着病情的发展和治疗手段的引入，纹理变化会更加明显，甚至出现变化扩散到其它的区域的现象。为了更准确的通过虹膜反映人体病变情况，就要针对虹膜的纹理变化进行自适应分区。

从虹膜解剖学我们知道，正常的卷缩论范围在虹膜图像的1/3。虹膜诊断学理论认为，卷缩轮范围变化与消化系统疾病密切相关。如卷缩轮紧缩表示胃部疾病，卷缩轮膨胀，反映着肠道系统疾病。又如：人眼正常虹膜不会出现裂缝，在慢性疾病病理中，会出现这种长长的裂缝。本发明提出的基于基元的串行区域生长与合并的虹膜图像动态分区技术，较好地解决了病理情况下虹膜病变特征的定位与测量问题。

串行区域生长与合并技术是分析图像中区域结构和相似性的模式识别方法，是建立在集合代数基础上，用集合论与概率论方法定量描述图像特征的技术。根据分割过程中处理策略的不同，分割技术又可分为并行技术与串行技术，虽然串行技术所需要计算时间比并行技术略长，但其抗噪声能力通常较强。

我们所提出的动态自适应虹膜分区(图谱覆盖)方法针对虹膜图像纹理复杂、区域变化不明显特性，在基于串行区域生长与合并技术基础上，通过对特征先验知识，定义区域生长与合并准则，以“基元”为区域生长点，这样在理论上满足医学图像处理的统计学要求，又可以便于计算机处理，通过实验得到了很好的结果。具体地方案实施如图5所示，包含如下步骤：

1.虹膜图像的采集；

2.虹膜图像的预处理；

3.虹膜定位；

4.CADIC虹膜图谱覆盖；

5.基于基元定义和区域生长的阳光放射沟提取；

6.虹膜自适应分区。

各个处理步骤的细节如下：

(1)虹膜图像的采集

虹膜具有对光的生理反应的特性；毋须物理接触，就可以在一定的距离获得被试的纹理图样。采样方法是对视频图像抓拍；采样者头部固定并正视；采样时所采用的光源产生的光斑打在瞳孔和虹膜的交界处；采样图像为768*576的24位位图。每个采样者采3到5组虹膜图像。

(2)虹膜图像的预处理

由于虹膜图像的采集受到光照及仪器差异的影响，包含了大量的噪声，同时眼睑和眼睫毛可能会掩盖一部分的虹膜，而睫状体对于虹膜诊病也是无用的信息，所以在对虹膜图像进行分区分析之前，要进行预处理工作。虹膜图像的预处理主要包括彩色图像灰度转换、噪声处理、几何校正、图像增强。

(2.1)彩色图像灰度化

结合图6-图7，目前，因为我们只对虹膜图像的纹理特征感兴趣，为了减少计算量，我们对得到的彩色图像进行灰度化。对采集到的彩色图像，应用如下公式将其变换为灰度图像。

$\left[\begin{array}{c}Y\\ I\\ Q\end{array}\right]=\left[\begin{array}{ccc}0.299& 0.587& 0.114\\ 0.596& -0.274& -0.322\\ 0.211& -0.523& 0.312\end{array}\right]\left[\begin{array}{c}R\\ G\\ B\end{array}\right]---(2-1)$

Y＝0.299×R+0.587×G+0.114×B            (2-2)

式中R——真彩图中每个像素的红分量；

G——真彩图中每个像素的绿分量；

B——真彩图中每个像素的蓝分量；

Y——该像素所对应的灰度

(2.2)虹膜图像灰度归一化

由于受到光照、拍照角度等因素的影响，使得图像中虹膜的尺寸和旋转的角度产生较大差异，这将极大的降低模板匹配的精度。虹膜归一化的目的是将原始虹膜图像中感兴趣的区域映射到相同的尺寸和相同的对应区域中，从而消除平移、旋转对纹理检测的影响。

(2.3)光斑填充

使用带有双光源的CCD采集设备采集到的虹膜图像，来自光源的入射光在舌面形成镜面反射，引起虹膜表面颜色的突变而形成光斑，且光斑颜色与光源颜色一致。因而虹膜图像上属于反光点的像素与虹膜上其他像素相比具有最高的亮度值和最低的饱和度值。这些光斑区域这会影响到定位结果，需要进行光斑填充。

结合图8-9，采用数字图像的自适应插值法来对光斑进行填充。零阶插值算法，即令输出像素的灰度值等于离它所映射的位置最近的输入像素的灰度值。

(3)虹膜定位

虹膜定位的目的就是除去那些于虹膜无关的信息(如眼睑Eyelid，睫毛Eyelash，瞳孔Pupil等)，也就是从虹膜图像中定位出虹膜的内边界和外边界。准确的虹膜定位是实现有效特征提取的前提。

(3.1)虹膜粗定位

结合图10，粗定位的目的是估算出虹膜内、外边缘对应圆的圆心和半径。

通过对虹膜图像的分析，得出以下的先验知识：瞳孔的直径最小也要大于120个像素；瞳孔是眼睛图像中面积最大的低灰度区域，而其圆心位置应该在上下左右灰度和值最小的区域附近，利用这些先验知识，可以减小搜索的空间提高精度。本文采用传统的灰度投影法对虹膜内、外边缘粗定位，得到虹膜内径r、外径r1和瞳孔中心(x_c，y_c)；

(3.2)虹膜精定位

结合图11，本发明采用Daugman算法来对虹膜进行精定位，其核心在于利用一个有效的积分微分算子来计算圆参数。虹膜具有良好的环状特性，所以可采用圆形如下算子对虹膜边界进行精定位。

在虹膜粗定位瞳孔中心(x_c，y_c)和半径r的基础上，采用圆探测算子实现虹膜精定位；所述虹膜精定位方法为基于圆探测算子的定位方法：

a、对于中心矩阵内的每一个元素，分别计算半径为r同一圆周上的像素点灰度和sum，其中中心矩阵的范围是(x_c-5，y_c-5)，半径变化范围为r-5～r+5；

b、对于r的变化值d，计算出对应的高斯函数值：

$g\left(d\right)=(1/(0.1\&times;\left(\sqrt{2\mathrm{\&pi;}}\right)))e^{((-{(d-6)}^2)/(2\&times;{0.1}^2))};$

c、计算灰度梯度值向量：det＝|g(d+1)-g(d)|×sum；

d、依次计算所有元素对应的det向量值；

e、在所有的det向量中取元素最大值，即为精定位的结果。

(4)CADIC虹膜图谱覆盖

结合图12，在虹膜内外边缘精确定位后，将虹膜部分所在区域的直角坐标(x，y)形式转换为极坐标(r，θ)形式。根据虹膜图谱和r、θ的值来将虹膜部分分成若干个区域。这是第一次分区，也就是理论分区。

(5)基于基元定义和区域生长的阳光放射沟提取

结合图13，分析虹膜主要是分析其结构纹理特征，虹膜诊病更是要将病灶特征独立而又完整地提取出来。本发明选取虹膜图像中疾病形成的阳光放射沟为处理对象，首先定义纹理基元并讨论基于基元的纹理表达方法。

设：纹理基元为h(x，y)，排列规则为r(x，y)，则纹理t(x，y)定义为

$t(x,y)=h(x,y)\&CircleTimes;r(x,y)---(4-1)$

设x_m，y_m为脉冲函数的位置坐标，则

r(x，y)＝∑δ(x-x_m，y-y_m)            (4-2)

根据卷积定理，在频域有：

T(u，v)＝H(u，v)R(u，v)              (4-3)

R(u，v)＝T(u，v)H(u，v)^-1            (4-4)

那么对于给定纹理基元h(x，y)的描述，即可得到反卷积滤波器H(u，v)^-1。将这个滤波器用于待处理纹理图像，每个脉冲都在纹理基元的中心。对整幅图像中不同纹理基元的分布进行统计可获得图像的全面纹理信息。

通过纹理基元确定对了对待分割的区域找到一个种子子图---即基元脉冲中心，然后将种子像素周围邻域中与种子像素有相同或相似性质的像素(通过定义的相同或相似规则)合并到种子像素所在区域中。在将新像素作为种子继续扩展，直到再没有满足条件的像素可被包括，这样的一个区域就长成。

对于一个给定的目标图像X和一个结构基元S，以符合S的种子区域

为中心在图像上进行迭代生长。

对于种子的每一个邻近像素x，在生长过程中应始终保持：

if：S(x)＝TRUE                (4-5)

$\mathrm{then}:X\&CirclePlus;S\&LeftArrow;X\&CirclePlus;S\&cup;\left\{x\right\}---(4-6)$

即，如果区域邻近像素x能符合结构基元定义，则将像素x合并到特征区域，使之成为新种子的一部分。

如果对于种子

的所有邻近像素x，均有：

S(x)＝FALSE                   (4-7)

则生长停止。

对于阳光放射沟特征，通过前期对临床病理的研究分析，从获得的先验统计知识定义像素生长规则如下。

规则1：沟长度为L个像素，5≤L≤R_虹膜；

规则2：沟宽度为W个像素，3≤W≤5；

规则3：沟的延长线与瞳孔交点的切线相垂直；

规则4：根据虹膜图像灰度直方图确定沟区域灰度最优阈值T_G；

规则5：沟的起点为卷缩轮边界点。

结合图14，根据CADIC虹膜分区图，对每个感兴趣的诊断候选区域，执行如下阳光放射沟检测算法。

1)计算阳光放射沟灰度阈值：

对整个虹膜统计灰度直方图，定义灰度阈值T_G，使得灰度小于T_G的像素数目为整个虹膜像素数目的10％。根据不同的虹膜图像采集设备，此比例数应根据人工标注阳光放射沟情况的先验统计数据做出相应的调整。

2)对任一个候选区Q，执行以下3)～10)。

3)沟起点选取：

对于任意像素x∈Q，如果有符合规则5的像素x的灰度值f(x)≤T_G，则认为x是一个沟起点。把找到的所有沟起点记为d_n，n＝1，2，3...N。

4)如果N＝0，即没有找到任何沟起点，则说明没有阳光放射沟存在，转10)结束本算法。

5)沟起点筛选与合并，形成种子候选：

对沟底点记d_n，n＝1，2，3...N，先后执行生长与合并操作，合并成种子候选d’_n，n＝1，2，3...N’，N’≤N。

6)种子候选筛选，形成种子：

对每个种子候选d’_n定义它在平行于虹膜内边缘方向上的像素跨度数为其宽度W(d’_n)，删除所有w(d’_n)＞5的种子候选d’_n。

所有剩余的种子候选d’_n成为种子，记为z_m，m＝1，2，3...M。在不造成混淆的情况下，也用z_m代表本种子的中心像素。

7)如果M＝0，即没有找到任何种子，则说明没有阳光放射沟存在，转10)结束本算法。

8)区域生长，形成阳光放射沟候选：

在每个种子z_m，m＝1，2，3...M周围邻域中搜索同种子像素具有相同或相似性质的像素，合并到种子像素所在区域中；在将新像素作为种子继续迭代扩展，直到再没有满足条件的像素可被包括，这样一个区域就长成了。

(8.1)定义每个种子(中心)z_m的邻域为NB(z_m)，NB(z_m)应符合规则1，2和3；

(8.2)定义基元结构为S，此处的S定义为典型阳光放射沟的横截线上的像素集合。

(8.3)对于任意像素x∈NB(z_m)，根据公式(4-5)、(4-6)所定义的生长规则进行区域生长。

9)阳光放射沟筛选：

对获得的阳光放射沟候选，根据规则1和2进行筛选。

10)输出结果，结束算法。

通过上述方法实验，可以获得满足裂缝基元定义的线条，即为阳光放射沟。

(6)虹膜动态自适应分区

具体做法体现在对病理特征(此处为阳光放射沟，也可以推广到其他虹膜特征上)在分区边界上的处理：当分区时病理特征被分割开来，不能被完整地表达在同一个区域中时，进一步从保持病理特征完整性的角度，重新确定该特征所覆盖区域。

6.1)首先，通过估计病理特征的区域进行病理特征筛选：

对病理特征边缘形状的提取采用估算的方法。求出组成特征的像素点在极坐标中的角度和极径的最大和最小值，这样就大概确定出了纹理所覆盖的区域的面积。从前面特征提取的结果可以看到，只有一部分重要的纹理特征(阳光放射沟)具有诊病价值。我们根据所确定的每个特征区域的面积，选取最重要的若干个特征，比如5个；删除其余的不重要特征。

6.2)结合图15，进行CADIC图谱覆盖，方法如前述。

6.3)判断是否有跨区域的病理特征：

对每一个虹膜图谱上的分区注意考察，如果发现其边界同重要特征有交叠，则记录这个交叠。

6.4)进行区域动态调整，完成自适应分区：

如果一个纹理特征覆盖多个区域，进行如下处理：对每个虹膜图谱分区，统计区域外内特征的像素点数目M；之后对其每处交叠边界，考察此区域外该特征的像素点并进行统计数目M’。

如果区域外像素点中只有极少部分M’＜＜M(例如M’小于5％的M)，则对这些少数像素点所分布的区域不予考虑；

如果区域外像素点过于重要M’＞M(说明该特征不是主要处于本区域)，则同样不予考虑，并删除本区域中的特征；

此外，说明区域外特征比较重要；则将此边界邻接的区域同本区域合并，实现了重新的区域划分，也就是进行了自适应分区。

例如在图16中跨区纹理特征进行自适应分区：在(l-1-2)、(l-5-2、l-5-3)与(l-7-2、l-7-3、l-7-4)区域中出现长的连续纹理特征---阳光放射沟，(l-1-2)区中特征没有出现跨区域现象；而在5区与7区中，把两个连续特征分割在几个小的区域中，使诊断很困难。应用上述发明技术，使图谱在覆盖虹膜同时，对检测到的完整特征实现自适应分区，确保特征信息的完整性，有利于诊病判断。

Claims (8)

1.一种基于基元结构定义与区域生长技术的虹膜分区及阳光放射沟提取方法，其特征在于在它是在基于串行区域生长与合并技术基础上，通过对特征先验知识，定义区域生长与合并准则，以基元为区域生长点，具体的方案包含如下步骤：(1)虹膜图像的采集；(2)虹膜图像的预处理；(3)虹膜定位；(4)计算机辅助诊断虹膜图谱CADIC虹膜图谱覆盖；(5)基于基元定义和区域生长的阳光放射沟提取；(6)虹膜自适应分区；

其中，

所述的基于基元定义和区域生长的阳光放射沟提取方法是选取虹膜图像中疾病形成的阳光放射沟为处理对象，首先定义纹理基元并讨论基于基元的纹理表达方法：

设：纹理基元为h(x，y)，排列规则为r(x，y)，则纹理t(x，y)定义为

$t(x,y)=h(x,y)\&CircleTimes;r(x,y)$

设x_m，y_m为脉冲函数的位置坐标，则

r(x，y)＝∑δ(x-x_m，y-y_m)

根据卷积定理，在频域有：

T(u，v)＝H(u，v)R(u，v)

R(u，v)＝T(u，v)H(u，v)^-1

对于给定纹理基元h(x，y)的描述，即能够得到反卷积滤波器H(u，v)^-1，将这个滤波器用于待处理纹理图像，每个脉冲都在纹理基元的中心，对整幅图像中不同纹理基元的分布进行统计能够获得图像的全面纹理信息；

通过纹理基元确定了对待分割的区域找到一个种子子图，即基元脉冲中心，然后将种子像素周围邻域中与种子像素有相同或相似性质的像素合并到种子像素所在区域中；再将新像素作为种子继续扩展，直到再没有满足条件的像素能够被包括，这样的一个区域就长成；

对于一个给定的目标图像X和一个结构基元S，以符合S的种子区域

为中心在图像上进行迭代生长；

对于种子

的每一个邻近像素x，在生长过程中应始终保持：

if：S(x)＝TRUE

then： $X\&CirclePlus;S\&LeftArrow;X\&CirclePlus;S\&cup;\left\{x\right\}$

即，如果区域邻近像素x能符合结构基元定义，则将像素x合并到特征区域，使之成为新种子的一部分；

如果对于种子

的所有邻近像素x，均有：

S(x)＝FALSE

则生长停止；

对于阳光放射沟特征，通过前期对临床病理的研究分析，从获得的先验统计知识定义像素生长规则如下：

规则1：沟长度为L个像素，5≤L≤R_虹膜；

规则2：沟宽度为W个像素，3≤W≤5；

规则3：沟的延长线与瞳孔交点的切线相垂直；

规则4：根据虹膜图像灰度直方图确定沟区域灰度最优阈值T_G；

规则5：沟的起点为卷缩轮边界点。

2.根据权利要求1所述的基于基元结构定义与区域生长技术的虹膜分区及阳光放射沟提取方法，其特征在于：虹膜图谱定义为：对于一幅虹膜图像，经过前期预处理与虹膜内外边缘精确定位后，将虹膜图像置于极坐标系(θ，r)中，把虹膜分为以瞳孔圆心为极点、以瞳孔边缘为起始边界的五个环、12个区的虹膜分区图谱，设一幅虹膜图像用A表示，R_T为瞳孔半径，R_I为虹膜半径，R_D为图像中卷缩轮的理论半径，R′_D为卷缩轮的实际半径，极坐标的极轴方向取为时钟的零点方向，极角按照顺时针方向定义，将A划分为五个环12个区后，其中的任意一个区用

来表示，其中：k＝1，2，...，12；i＝1，2，...，5，每个分区的定义如下：

$A_i^k(\mathrm{\&theta;},r):\frac{(k-1)\mathrm{\&pi;}}{6}<\mathrm{\&theta;}\&le;\frac{\mathrm{k\&pi;}}{6}$

$\left\{\begin{array}{c}\mathrm{if}:i=1;\mathrm{then}:r=R_T\\ \mathrm{if}:i=2;\mathrm{then}:R_T<r\&le;R_T+\frac{2}{5}(R_I-R_T)\\ \mathrm{if}:i=\mathrm{3,4,5};\mathrm{then}:R_T+(i-1)\frac{(R_I-R_T)}{5}<r\&le;R_T+i\frac{(R_I-R_T)}{5}\end{array}\right..$

3.根据权利要求1所述的基于基元结构定义与区域生长技术的虹膜分区及阳光放射沟提取方法，其特征在于：虹膜图像的采集方法是对视频图像抓拍，采样者头部固定并正视，采样时所采用的光源产生的光斑打在瞳孔和虹膜的交界处，采样图像为768*576的24位位图，每个采样者采3到5组虹膜图像。

4.根据权利要求3所述的基于基元结构定义与区域生长技术的虹膜分区及阳光放射沟提取方法，其特征在于：虹膜图像的预处理方法是在对虹膜图像进行分区分析之前进行预处理工作，虹膜图像的预处理方法包括彩色图像灰度转换、噪声处理、几何校正、图像增强方法，图像灰度转换方法为对采集到的彩色图像，应用如下公式将其变换为灰度图像：

$\left[\begin{array}{c}Y\\ I\\ Q\end{array}\right]=\left[\begin{array}{ccc}0.299& 0.587& 0.114\\ 0.596& -0.274& -0.322\\ 0.211& -0.523& 0.312\end{array}\right]\left[\begin{array}{c}R\\ G\\ B\end{array}\right]$

Y＝0.299×R+0.587×G+0.114×B

式中R——真彩图中每个像素的红分量；

G——真彩图中每个像素的绿分量；

B——真彩图中每个像素的蓝分量；

Y——该像素所对应的灰度；

再进行虹膜图像灰度归—化，通过将原始虹膜图像中感兴趣的区域映射到相同的尺寸和相同的对应区域中进行虹膜归一化，然后进行光斑填充，采用数字图像的自适应插值法来对光斑进行填充，零阶插值算法，即令输出像素的灰度值等于离它所映射的位置最近的输入像素的灰度值。

5.根据权利要求1所述的基于基元结构定义与区域生长技术的虹膜分区及阳光放射沟提取方法，其特征在于虹膜定位方法是除去那些于虹膜无关的信息，也就是从虹膜图像中定位出虹膜的内边界和外边界，包括：

(1)虹膜粗定位：采用灰度投影法对虹膜内、外边缘粗定位，得到虹膜内径r、外径r1和瞳孔中心(x_c，y_c)；

(2)虹膜精定位：采用Daugman算法来对虹膜进行精定位，利用一个有效的积分微分算子来计算圆参数，采用圆形探测算子对虹膜边界进行精定位，在虹膜粗定位瞳孔中心(x_c，y_c)和半径r的基础上，采用圆探测算子实现虹膜精定位；所述虹膜精定位方法为基于圆探测算子的定位方法：a、对于中心矩阵内的每一个元素，分别计算半径为r同一圆周上的像素点灰度和sum，其中，中心矩阵的范围是(x_c-5，y_c-5)，半径变化范围为r-5～r+5；b、对于r的变化值d，计算出对应的高斯函数值： $g\left(d\right)=(1/(0.1\&times;\left(\sqrt{2\mathrm{\&pi;}}\right)))e^{((-{(d-6)}^2)/(2\&times;{0.1}^2))};$ c、计算灰度梯度值向量：det＝|g(d+1)-g(d)|×sum；d、依次计算所有元素对应的det向量值；e、在所有的det向量中取元素最大值，即为精定位的结果。

6.根据权利要求1所述的基于基元结构定义与区域生长技术的虹膜分区及阳光放射沟提取方法，其特征在于：计算机辅助诊断虹膜图谱CADIC虹膜图谱覆盖的方法是在虹膜内外边缘精确定位后，将虹膜部分所在区域的直角坐标(x，y)形式转换为极坐标(r，θ)形式，根据虹膜图谱和r、θ的值来将虹膜部分分成若干个区域。

7.根据权利要求1所述的基于基元结构定义与区域生长技术的虹膜分区及阳光放射沟提取方法，其特征在于：根据计算机辅助诊断虹膜图谱CADIC虹膜分区图，对每个感兴趣的诊断候选区域，执行如下阳光放射沟检测算法：

1)计算阳光放射沟灰度阈值：

对整个虹膜统计灰度直方图，定义灰度阈值T_G，使得灰度小于T_G的像素数目为整个虹膜像素数目的10％，根据不同的虹膜图像采集设备，此比例数应根据人工标注阳光放射沟情况的先验统计数据做出相应的调整；

2)对任一个候选区Q，执行以下3)～10)；

3)沟起点选取：

对于任意像素x∈Q，如果有符合规则5的像素x的灰度值f(x)≤T_G，则认为x是一个沟起点，把找到的所有沟起点记为d_n，n＝1，2，3...N；

4)如果N＝0，即没有找到任何沟起点，则说明没有阳光放射沟存在，转10)结束本算法；

5)沟起点筛选与合并，形成种子候选：

对沟底点记d_n，n＝1，2，3...N，先后执行生长与合并操作，合并成种子候选d’_n，n＝1，2，3...N’，N’≤N；

6)种子候选筛选，形成种子：

对每个种子候选d’_n定义它在平行于虹膜内边缘方向上的像素跨度数为其宽度W(d’_n)，删除所有W(d’_n)＞5的种子候选d’_n；

所有剩余的种子候选d’_n成为种子，记为z_m，m＝1，2，3...M，在不造成混淆的情况下，也用z_m代表本种子的中心像素；

7)如果M＝0，即没有找到任何种子，则说明没有阳光放射沟存在，转10)结束本算法；

8)区域生长，形成阳光放射沟候选：

在每个种子z_m，m＝1，2，3...M周围邻域中搜索同种子像素具有相同或相似性质的像素，合并到种子像素所在区域中；在将新像素作为种子继续迭代扩展，直到再没有满足条件的像素能够被包括，这样一个区域就长成了；

(8.1)定义每个种子z_m的邻域为NB(z_m)，NB(z_m)应符合规则1，2和3；

(8.2)定义基元结构为S，此处的S定义为典型阳光放射沟的横截线上的像素集合；

(8.3)对于任意像素x∈NB(z_m)，根据上述定义的生长规则进行区域生长；

9)阳光放射沟筛选：对获得的阳光放射沟候选，根据规则1和2进行筛选；

10)输出结果，结束算法；

通过上述方法实验，能够获得满足裂缝基元定义的线条，即为阳光放射沟。

8.根据权利要求7所述的基于基元结构定义与区域生长技术的虹膜分区及阳光放射沟提取方法，其特征在于：虹膜动态自适应分区方法是在对病理特征在分区边界上的处理：当分区时病理特征被分割开来，不能被完整地表达在同一个区域中时，进一步从保持病理特征完整性的角度，重新确定该特征所覆盖区域：

1)首先，通过估计病理特征的区域进行病理特征筛选：

对病理特征边缘形状的提取采用估算的方法，求出组成特征的像素点在极坐标中的角度和极径的最大和最小值，这样就大概确定出了纹理所覆盖的区域的面积，从前面特征提取的结果能够看到，只有一部分重要的纹理特征具有诊病价值，根据所确定的每个特征区域的面积，选取最重要特征，删除其余的不重要特征；

2)进行计算机辅助诊断虹膜图谱CADIC图谱覆盖，方法如前述；

3)判断是否有跨区域的病理特征：

对每一个虹膜图谱上的分区注意考察，如果发现其边界同重要特征有交叠，则记录这个交叠；

4)进行区域动态调整，完成自适应分区：如果一个纹理特征覆盖多个区域，进行如下处理：对每个虹膜图谱分区，统计区域内特征的像素点数目M；之后对其每处交叠边界，考察此区域外该特征的像素点并进行统计数目M；

如果区域外像素点中只有极少部分M’＜＜M，M’小于5％的M，则对这些少数像素点所分布的区域不予考虑；

如果区域外像素点过于重要M’＞M，则同样不予考虑，并删除本区域中的特征；

此外，说明区域外特征比较重要；则将此边界邻接的区域同本区域合并，实现了重新的区域划分，也就是进行了自适应分区。

Images