一种基于图像配准的视网膜血管壁厚度测量方法与流程

本发明涉及图像处理技术领域，具体涉及一种基于图像配准的视网膜血管壁厚度测量方法。

背景技术

眼底视网膜血管是人体中唯一可看见的血管，医生把它当作了解其它脏器血管情况的窗口。如高血压病人眼底可见到视网膜动脉硬化、糖尿病病人眼底可见毛细血管瘤、小的出血点和渗出物等。眼底血管管壁厚度能作为糖尿病视网膜病变的早期指标以及为判断血管静脉硬化提供帮助，而眼底彩色图像包含血管壁外侧、眼底荧光造影图像仅包含血管内侧。配准眼底图像能够使两幅图的信息互补，并使得测量血管壁厚度更加容易。

图像配准是指将同一场景的不同图像“对齐”或进行广义的配准，以消除存在的几何畸变。对同一场景使用相同或不同的传感器，在不同条件下获取的两个或多个图像一般都会存在差异。同一场景的多次成像的差别可以表现在：不同的分辨率、不同的灰度属性、不同的位置(平移和旋转)、不同的尺度、不同的非线性变换等等。多种成像模式产生的图像(称之为多模态图像)会表现出不同的分辨率、不同的灰度属性等差异。图像配准的主要目的是消除几何上的差异，对于灰度信息的差异并不关心。

目前测量眼底血管壁厚度的方法为直接通过眼底彩色图像中血管的两侧来测量，但是，眼底彩色图像中血管的两侧是反光形成的，并不是真正的血管壁厚度，因此，采用这种直接测量视网膜血管壁厚度的方法导致方法测量结果不准确，影响医生对患者是否含有糖尿病等疾病的判断。

此外，视网膜血管壁的平均厚度，能够反应糖尿病、脑卒中等疾病发生情况，而单个血管壁厚度还能反应该血管是否出现坏死、血管壁破裂等症状。这是因为单个血管壁厚度出现问题的时候，没有影响到平均血管壁厚度而带来的误判。因此对于视网膜血管壁厚度，既要测量平均厚度、也要测量单个血管壁厚度。

技术实现要素：

本发明的目的在于克服现有技术的不足，提供一种基于图像配准的视网膜血管壁厚度测量方法，通过输入同一眼球的一张眼底彩色图像和一张眼底荧光造影图像，配准两幅图像，提高了配准的准确率；此外通过标注视盘后计算出眼底血管壁平均厚度及单个血管壁厚度，提高了血管壁厚度测量结果的准确性，有助于医生进行判断，节约了医生的时间。

为实现上述发明目的，本发明一种眼底图像配准的视网膜血管壁厚度测量方法，包括以下步骤：

s1：获取眼底彩色图像和眼底荧光造影图像；

s2：对获取的眼底彩色图像和眼底荧光造影图像进行预处理；

s3：利用pca算法对预处理后的眼底彩色图像和眼底荧光造影图像进行血管分割；

s4：利用图像腐蚀和扩张方法对分割后的眼底彩色图像和眼底荧光造影图像进行去噪处理；

s5：利用zhang并行快速细化算法对去噪处理后的眼底彩色图像和眼底荧光造影图像进行骨骼化处理；

s6：对骨骼化处理后的眼底彩色图像和眼底荧光造影图像进行harris角点检测，得到包含分支点和交叉点的角点；

s7：以harris角点检测得到的各角点为中心，判断骨骼化处理后的眼底彩色图像和眼底荧光造影图像中的各角点是分支点还是交叉点，若是交叉点，则舍弃该角点，若是分支点，则从该分支点的三个分支上分别找下一个角点，若三个分支上的三个角点同样是含有三个分支的角点，则把该角点作为特征点中心，将分支上的三个角点和该角点组合成特征点，特征点中心以外的角点采用顺时针记录；分别记录眼底彩色图像及眼底荧光造影图像中的特征点并构成两个特征点集；

s8：将眼底彩色图像的特征点集中的每个角点分别匹配眼底荧光造影图像的特征点集中的每个角点，将每个特征点中的角点横、纵坐标代入仿射变换公式并求得参数集，根据参数集中的参数计算匹配度，将参数集中匹配度最高的参数作为最终进行图像配准的仿射变换参数；

s9：以步骤s1获取的眼底荧光造影图像作为参考图像，以眼底彩色图像作为浮动图像，根据步骤s8得到的仿射变换参数，将获取的眼底彩色图像和眼底荧光造影图像进行配准；

s10：对配准后的图像进行视盘标注，并以视盘中心为圆心，视盘直径dd为基准，得到直径为2dd到3dd的圆构成的圆环；

s11：分别计数配准后图像中圆环内的像素点总数或某一血管段的像素点总数，根据像素点总数和视盘直径dd，计算视网膜血管壁厚度，所述视网膜血管壁厚度包括视网膜中单个血管的血管壁厚度p和视网膜血管壁的平均厚度p。

本发明的有益效果是：由步骤s8中提出的新的配准方法以及步骤s7中新的特征点提取方法提高了配准的准确率；步骤s9中对原始图像进行配准嵌合并生成的图像，能够辅助医生进行判断和提升医生判断病情的准确性及节约分析病情的时间。步骤s10、s11中新的血管壁厚度测量方法提高了眼底血管壁厚度的准确度，通过平均厚度来判断患者是否含有糖尿病等疾病，以及单个厚度判断血管是否有血管壁破裂等情况，可大大节约医生的时间，对眼科临床有较大益处。

附图说明

图1是本发明基于图像配准的视网膜血管壁厚度测量方法的流程图；

图2是本发明图像预处理结果示意图；

(a)为眼底彩色图像示意图；(b)为眼底荧光造影图像示意图；

图3是本发明血管分割结果示意图；

(a)为眼底彩色图像示意图；(b)为眼底荧光造影图像示意图；

图4是本发明图像骨骼化处理示意图；

(a)为眼底彩色图像示意图；(b)为眼底荧光造影图像示意图；

图5是本发明特征点匹配结果示意图；

图6是本发明图像配准结果示意图；

图7是本发明视盘标注示意图。

具体实施方式

下面结合附图对本发明的具体实施方式进行描述，以便本领域的技术人员更好地理解本发明。需要特别提醒注意的是，在以下的描述中，当已知功能和设计的详细描述也许会淡化本发明的主要内容时，这些描述在这里将被忽略。

为了解决目前直接测量视网膜血管壁厚度的方法存在的结果不准确，影响医生对疾病的判断的问题，本发明提成了一种基于多模态眼底图像配准的视网膜血管壁厚度测量方法，包括以下步骤：

s1：获取眼底彩色图像和眼底荧光造影图像；获取的眼底彩色图像与眼底荧光造影图像均来自同一个人的同一个眼球，所述眼底彩色图像与眼底荧光造影图像均包括图像宽、高与颜色分量信息，且图像宽、高一致。

s2：对眼底彩色图像和眼底荧光造影图像进行预处理，具体为：

本步骤中，如果是眼底彩色图像，则进行通道处理和/或灰度处理，这里如果图像模糊或信息不全时，则要进行通道处理和灰度处理，经过处理后图像准确度提高了，但会降低效率，若图像清晰，则进行通道处理或灰度处理，其中通道处理是提取rgb三个色域中的一个，灰度处理则是将图片处理成灰度图像，灰度值的范围为0-255，表示亮度从深到浅，颜色从黑到白。若是眼底荧光造影图像，则进行反色处理，反色处理则是将灰度图像的每个值减去255后求绝对值。

s3：利用pca算法对预处理后的眼底彩色图像和眼底荧光造影图像进行血管分割；

本步骤中，将眼底彩色图像和眼底荧光造影图像的rgb三通道中任意两通道分量置0，对非零通道分量向量化后进行pca分析，得到score矩阵；

对score矩阵中非零通道分量归一化，然后进行clahe直方图均衡；

对均衡结果进行中值平滑，然后与平滑前图像做差，将差值与阈值比较，得到二值化图像，利用bwareaopen函数去除二值化图像中面积小的区域，然后取反色并将整体背景设为黑色，得到分割图像。

s4：利用图像腐蚀和扩张方法对分割后的眼底彩色图像和眼底荧光造影图像进行去噪处理；

本步骤中，首先对图像进行腐蚀和膨胀处理操作，这样做可以连通血管间分割效果较差的分支，另外可以去除分割时把一条血管分割成两条的错误情况。此时血管中有部分点是无法连通，只要去掉这些无法连通的点即可完成去噪。

s5：利用zhang并行快速细化算法对去噪处理后的眼底彩色图像和眼底荧光造影图像进行骨骼化处理；

(1)遍历去噪处理后的眼底彩色图像和眼底荧光造影图像中的每个血管像素点，并以每个血管像素点pi为中心，构成窗口w：

其中，i为第i个血管像素点，i＝1、2……n-8，n为血管像素点总数，pi+1～pi+8为血管像素点pi周围的8个像素点，若pi+1～pi+8中的像素点为血管像素点，则值为1，否则值为0；令n(pi)表示pi周围的像素点值之和，s(pi)表示按pi～pi+8排列出现0-1的次数；

(2)第一次遍历所有的血管像素点pi，判断窗口w是否满足2<＝n(pi)<＝6、s(pi)＝1、pi+1*pi+3*pi+5＝0、pi+3*pi+5*pi+7＝0的条件，若满足条件则将pi置为0，且将该血管像素点置为背景点，即不再是血管像素点；

(3)再次遍历剩下的所有血管像素点，若满足2<＝n(pi)<＝6、s(pi)＝1、pi+1*pi+3*pi+7＝0、pi+1*pi+5*pi+7＝0的点，则置为背景点；

(4)返回步骤(2)直到没有血管像素点置为背景点为止，得到二值图像细化后的骨架。

s6：对骨骼化处理后的眼底彩色图像和眼底荧光造影图像进行harris角点检测，角点包括分支、交叉点；利用移动的水平窗口和竖直窗口在骨骼化后的图像中的每个血管点移动来计算灰度变化值，用上、下、左、右四个方向移动后的窗口灰度值减去某个血管像素点的窗口灰度值即得到灰度变化值，若窗口在水平和竖直方向上移动后灰度值变化均高于某一阈值，则认为在窗口内遇到了角点即该血管像素点是角点，反之若灰度变化值在阈值之内，则认为该点并非角点；

其水平窗口为：

竖直窗口为：

对于以某一血管点pi为中心，构成的窗口w为：

其中，左移即w+sp，右移即w-sp，上移即w-sz、下移即w+sz。

s7：以harris角点检测得到的每个角点为中心，判断骨骼化处理后的眼底彩色图像和眼底荧光造影图像中的该角点是分支点还是交叉点，若是交叉点，则舍弃该角点，若是分支点，则从该分支点的三个分支上分别找下一个角点，若三个分支上的三个角点同样是含有三个分支的角点，则把该角点作为特征点中心，将分支上的三个角点和该角点组合成特征点；除特征点中心以外的角点采用顺时针记录；分别记录眼底彩色图像及眼底应该造影图像中的特征点构成两个特征点集；

s8：将s7中眼底彩色图像的特征点集中的每个点分别匹配眼底荧光造影图像的特征点集中的每个点，如首先分别取出一个角点，眼底彩色图像的特征点集中的点作为浮动点，眼底荧光造影图像的特征点集中的点作为参考点。由于一个特征点包含4个角点，且包含顺时针记录的3个角点，因此每次需要做3次匹配，每次匹配包括浮动点中心角点对应参考点中心角点，浮动点的3个顺时针角点对应参考点的3个顺时针角点，而每次匹配就是将角点的横纵坐标x、y、x'、y'带入仿射变换公式并求得参数向量，因为每次匹配都能得出3个参数(两次旋转)，因此参数集个数为特征点集乘积的3倍；

本步骤中将每个特征点中的角点横、纵坐标代入仿射变换公式并求得参数集，其中，仿射变换公式为：

其中a1～a4分别为cos(θ)、-sin(θ)、sin(θ)、cos(θ)，其中θ是仿射变换的角度参数，a1、a2、a3、a4、tx和ty分别为仿射变换公式参数，x、y分别为眼底荧光造影图像特征点中角点的横、纵坐标，x'、y'分别为眼底彩色图像特征点中角点的横、纵坐标，x、y、x'、y'均为已知值；

根据参数集中的参数计算匹配度，将参数集中匹配度最高的参数作为最终进行图像配准的仿射变换参数，其中，根据参数集中的参数计算匹配度有两种方式：

第一种：设置特征点之间的差值阈值，计算眼底彩色图像中的特征点经过仿射变换公式变换后在阈值范围内的个数，最后保留使得特征点匹配的数量最大的仿射变换公式参数，从而确定仿射变换参数。

第二种：直接遍历参数，将眼底彩色图像血管图像中每个血管像素点经过仿射变换映射到眼底荧光造影图像血管图像中，若映射的该像素点同样是血管像素点，则计数增加一，最后保留使得血管像素点匹配的计数最大的仿射变换公式参数，从而确定仿射变换参数。

从计数和时间复杂度角度看，第二种方法效率不如第一种方法，但更加准确，能够应对病变严重的眼底图像配准(血管不明显)。但提升的效果需要效率权衡，具体视实验需要效率还是精准度而定。

s9：以步骤s1获取的眼底荧光造影图像作为参考图像，以眼底彩色图像作为浮动图像，根据步骤s8得到的仿射变换参数，将获取的眼底彩色图像和眼底荧光造影图像进行配准。

本步骤中，将仿射变换参数和浮动图像中的每个像素点代入仿射变换公式，若变换后的像素点超过参考图像的范围，则抛弃该像素点；若在范围内，则保存该像素点，得到一组像素点，将该组像素点映射到参考图像的像素点数组中，则完成配准。

s10：对配准后的图像进行视盘标注，并以视盘中心为圆心，视盘直径dd为基准，得到直径为2dd到3dd的圆构成的圆环；

s11：分别计数配准后图像中圆环内的像素点总数，根据像素点总数的差值和视盘直径dd计算视网膜血管壁厚度。

本发明中视网膜血管壁厚度为视网膜中某血管壁厚度或视网膜血管壁的平均厚度方法，其计算方法分别如下：

第一种是：视网膜中某血管壁厚度，其具体过程为:

分别在眼底荧光造影图像和眼底彩色图像的圆环内某一血管段中设置探针标记点；

以该探针标记点为中心，向邻近像素点扩散，每扩散一个像素点，计数值count1、count2就分别增加一，直到扩散的所有像素点均在圆环之外时结束；此时计数值count1、count2的值则为该眼底荧光造影图像和眼底彩色图像中该血管的像素点总数，根据像素点总数差值和视盘直径dd计算该血管的血管壁厚度p，其计算公式为：

第二种是：视网膜血管壁的平均厚度方法，其具体过程为:

计算眼底荧光造影图像和眼底彩色图像的圆环内的血管内像素点总数为count_total1及count_total2，根据像素点总数差值和视盘直径dd计算视网膜平均血管壁厚度为，其计算公式为：

实施例1

本实施例先作如下说明：

将眼底彩色图像中用到的名称以序号1来标注，将眼底荧光造影图像中用到的名称以序号2来标注；即1代表眼底彩色图像，2代表眼底荧光造影图像。

实施例中的数组均可保存为图像，另外参加下一步骤中计算的数组均可来源于上一步骤中的结果，也可以重新读取上一步骤中生成的结果图像。此外以下数组及图片信息保存结果均以a到z的顺序保存，而不具有实际意义。

如图1所示，本发明基于多模态眼底图像配准的视网膜血管壁厚度测量方法，包括以下步骤：

s101：分别输入眼底彩色图像及眼底荧光造影图像；

本步骤中输入的眼底彩色图像a1与眼底荧光造影图像a2必须来自同一个人的同一个眼球，可以为不同时期。另外最好保证输出的a1、a2的宽、高一致，即像素大小一致，如均为768*768。输入图像后分别读取a1、a2中的信息，包括图像宽、高与颜色分量，将信息分别保存到五维数组中arr1_a、arr2_a中：

arr1_a＝[x1，y1，r1，g1，b1]；arr2_a＝[x2，y2，r2，g2，b2]

其中xi表示横坐标、yi表示纵坐标、ri表示红色分量、gi表示绿色分量、bi表示蓝色分量，i＝1，2。

s102：对眼底彩色图像和眼底荧光造影图像进行预处理；

本步骤中，如果是眼底彩色图像a1，则进行通道处理和/或灰度处理，这里如果图像模糊或信息不全时，则要进行通道处理和灰度处理，经过处理后图像准确度提高了，但会降低效率，若图像清晰，则进行通道处理或灰度处理，其中通道处理是提取rgb三个色域中的一个，灰度处理则是将图片处理成灰度图像，灰度值的范围为0-255，表示亮度从深到浅，颜色从黑到白。若是眼底荧光造影图像a2，则进行反色处理，反色处理则是将灰度图像的每个值减去255后求绝对值。

将s101中得到的arr1_a数组进行红色通道提取：

arr1-b＝[x1，y1，r1，r1，r1]；

可见g、b分量的值被置为r分量的值。

对s101中得到的arr2_a数组进行反色处理：

arr2_b＝255-arr2_a[x2，y2]＝[x2，y2，255-r2，255-g2，255-b2]；

此处用255减去分量值，或由分量值减去255取绝对值均可。

得到的新的arr1_b与arr2_b可保存为图像b1、b2作为结果展示，结果如图2所示。

s103：利用pca算法对预处理后的两幅血管图像进行血管分割；

本步骤中，血管分割包括：

(1)读取b1、b2图像生成的arr1_b与arr2_b数组，将其中的r、g、b任意两个分量值置为0；此时，将五维数组降维为三维数组arr1_c[x1][y1][l1]、arr2_c[x2][y2][l2]，其中l1、l2指r、g、b分量中非0的分量：

若保留红色分量，则三维数组值为：

arr1_c＝[x1，y1，r1]＝[x1，y1，l1]；

arr2_c＝[x2，y2，255-r2]＝[x2，y2，l2]；

若保留绿色分量，则三维数组为：

arr1_c＝[x1，y1，g1]＝[x1，y1，l1]；

arr2_c＝[x2，y2，255-g2]＝[x2，y2，l2]；

(2)分别对l1、l2分量向量化后进行pca分析：

特征中心化，即每一维的数据都减去该维的均值由于图片信息是分横坐标x，纵坐标y，以及分量l放入数组中的，如数组包含512行和768列，即x＝512、y＝768，那么共有x*y个像素点，均值就是这x*y个像素点的分量总和除以x*y，这里分别代入arr1_c、arr2_c数组：

其中，x1、y1为arr1_c数组的二维长度；x2、y2为arr2_c数组的二维长度。若y1、y2分辨率维512*768，此时x1、x2的值为512，y1、y2的值为768。i1、i2分别为l1、l2中的值。

(3)计算协方差矩阵c1、c2：

其中y1、y2为上述步骤中的值，li、1j为属于x1的值，表示第li、1j行，2i、2j为属于x2的值，表示第2i、2j行，相当于对arr1_cn、arr2_cn数组的每一行做协方差计算，如arr1_cn、arr2_cn装填的是512*768的图片，那么每一行都有768列，对于arr1_cn，分别取第0行和第1行做协方差运算时，c1中的li就是0，1j就是1。

分别计算协方差矩阵c1、c2的特征值和特征向量，并分别选取大的特征值对应的特征向量，得到新的数据集d1、d2；

(4)对数据集d1、d2分别进行归一化，然后进行clahe直方图均衡：图像分块，以块为单位，先计算直方图，然后修剪直方图，最后均衡；

(5)对均衡结果进行中值平滑，然后与平滑前图像做差，得到细节图像z1、z2；然后对细节图像z1、z2进行阈值处理，得到二值化图像bw1、bw2，利用bwareaopen函数分别将二值化图像bw1、bw2中面积小的区域去除，得到图像bc1、bc2，然后再取反色并将整体背景设为黑色，即可得到眼底图像血管分割后的图像f1、f2。

说明：上述步骤中的d1、d2、z1、z2、bw1、bw2、bc1、bc2、f1、f2均为图像数组，可将图像数组保存为图片，若不保存为图像，数组的处理结果分别为arr1_f、arr2_f：

此外f1、f2的图像保存结果如图3所示。

s104：利用图像腐蚀和膨胀方法对两幅血管图像进行去噪处理；

首先对f1、f2图像进行腐蚀和膨胀处理操作，这样做可以连通血管间分割效果较差的分支，另外可以去除分割时把一条血管分割成两条的错误情况。此时血管中有部分点是无法连通，只要去掉这些无法连通的点即可完成去噪，具体实施为：

代入s103血管分割后或读取f1、f2图片生成的数组arr1_f、arr2_f，调用腐蚀、膨胀函数，得到新的血管分割数组arr1_fn、arr2_fn并保存为图像fn1、fn2。

s105：利用zhang并行快速细化算法对fn1、fn2两幅血管图像进行血管的骨骼化操作，得到两幅血管骨架图像；由于一次提取会导致骨骼化不够平滑(非线性)，因此需要迭代提取骨架，具体实施步骤如下：

(1)分别遍历s104结果数组或读取fn1、fn2图片生成的arr1_fn、arr2_fn中所有前景像素点，此处前景像素点即为血管点，背景像素点为非血管像素点。

对于某个前景像素点pi，其周围8个像素点构成的窗口w为：

对符合待删除条件的像素点标记为删除；其中i为当前像素点，pi+1～pi+8为pi相邻的8个像素点，条件如下：

其中n(pi)表示和pi相邻的8个像素点中为前景像素点的个数；s(pi)表示从pi+1～pi+8像素点中出现0-1的累计次数(如pi为0、pi+1为1则构成一个0-1，pi+2为0、pi+3为0则不构成0-1)，其中0表示背景，1表示前景；

(2)与(1)中条件类似，满足条件的像素点p1则标记为删除，条件如下：

(3)循环上述两步骤，直到两步中都没有像素点被标记为删除为止，输出的结果即为二值图像细化后的骨架。

(4)按顺序将步骤(1)至步骤(3)封装为skeleton函数，分别代入arr1_fn、arr2_fn，迭代多次skeleton函数，直到间隔调用skeleton函数后数组没有变化，此时分别保存骨骼化数组为arr1_g、arr2_g。

(5)保存arr1_g、arr2_g数组为骨骼化图像g1、g2，如图4所示。

s106：对骨骼化后的两幅血管骨架图像g1、g2进行harris角点检测，得到包含分支点和交叉点的角点；

本步骤中，利用移动的窗口在图像中计算灰度变化值，其流程包括转化为灰度图像、计算差分图像、高斯平滑、计算局部极值、确认角点。对第一次的角点集合再进行候选点剔除，这是由于第一次检测出的角点不一定全是正确的分支点或交叉点，具体实施步骤如下：

(1)设定3x3的水平方向移动窗口：

和垂直方向移动窗口：

若调用水平方向移动窗口做处理，则为某一像素点的左右移动；若调用垂直方向移动窗口则为上下移动。若同时使用，则能形成左上、左下、右上、右下的移动。

(2)取窗口函数w(x,y)为矩形窗rectangle函数，其中，在窗口内的值为1，不在窗口内的值为0：

(3)从初始位置x＝0，y＝0开始移动窗口，其中平移量为[u,v]，得出灰度变化值；灰度变化公式为：

其中i(x+u,y+v)表示窗口移动后的灰度值，i(x,y)为某一血管像素点及周围8个像素点构成的窗口的灰度值。分别代入s105中结果数组或由g1、g2生成的arr1_g、arr2_g数组，遍历数组中的血管像素点，对于某一点及其周围的八个点构成的窗口w，左移即w+sp，右移即w-sp，上移即w-sz、下移即w+sz。

分别计算这上、下、左、右四种移动后的灰度变化值，即用现在窗口灰度值减去窗口w的灰度值取绝对值。若四个方向移动的灰度值变化均高于某个阈值，则认为该点为角点。

(4)按顺序将步骤(1)至步骤(3)封装为findharris检测函数，代入s105的结果数组arr1_g、arr2_g，输出被判定为角点的数组arr1_h、arr2_h：

s107：遍历arr1_h、arr2_h中每个角点，从骨骼化处理后的眼底彩色图像和眼底荧光造影图像中提取特征点，并匹配特征，得到仿射变换参数；具体实施步骤如下：

(1)输入s105的骨骼化数组arr1_g、arr2_g和s106得到的被盘点为角点的数组arr1_h、arr2_h，分别遍历arr1_h、arr2_h数组，从骨骼化数组arr1_g、arr2_g找到仅有三个分支的点f1，由于少于三个分支的点可能是误分的分支分叉点，而多余三个分支的点则认为是交叉点，这里暂不使用交叉点。

(2)沿着f1点的三个分支，查找下一个角点，由于有三个分支，因此可以找到三个分支点f2～f4，若这三个分支无法找到三个角点，则认为f1不是候选点；若分别查找到三个角点，则查看f2～f4是否为同样为分支分叉点，只有当它们同样是分支分叉点时才保存f1～f4到特征向量中，其中特征向量为s：

s＝[f1，f2，f3，f4，l1，l2，l3]

其中f1～f4包含f1～f4四个点的角度信息，由于都含有三个角度，因此：

f1＝[a1，a2，a3]

f2～f4以内类推，另外11～l3表示点f1分别到f2、f3、f4的距离。此外f2～f4构成顺时针结构。

(3)分别保存眼底荧光造影图像及眼底彩色图像的骨骼化图像中的所有特征点的特征向量s1、s2到数组fa1[s1]、fa2[s2]中，遍历fa1数组，以每个特征点来遍历fa2数组中的每个特征点，由于每个特征包含四个点，将这四个点分别代入仿射变换公式：

其中a1～a4分别为cos(theta)、-sin(theta)、sin(theta)、cos(theta)，其中theta是仿射变换的角度参数，此处确定theta的角度即可得出a1～a4的值，如theta＝60°，则a1＝1/2；x、y为参考图像特征点中角点的横纵坐标，x’、y’表示浮动图像特征点中角点的横纵坐标，此步骤旨在求取a1～a4、tx、ty的参数值，x、y、x’、y’为已知值。

这样求出a1、a2、a3、a4、tx、ty的值，即得到仿射变换公式参数；若遍历所有的特征向量，可得到fa1数组长度乘以fa2数组长度组这样的仿射变换公式，将仿射变换公式参数集保存到数组fang中：

fang＝[a1，a2，a3，a4，tx，ty]；

(4)确定最终仿射变换参数

这里有两种方式：

第一种：设置特征点之间的差值阈值如20，即一个特征点中的四个角点与另一个特征点中的四个角点之间的距离最小值之和，遍历数组fang，求取眼底彩色图像中的特征点经过仿射变换后在阈值范围内的个数，数量最多的仿射变换公式参数则为最后保留的参数fang_final，特征点匹配结果如图5所示。

第二种：直接遍历数组fang，将眼底彩色图像的血管分割图像的每个血管像素点经过仿射变换映射到眼底荧光造影图像的血管分割图像中，若映射的该点同样是血管像素点，则计数增加一。最后同样保留使得计数最大的仿射变换公式参数fang_final。

两种方法同样可行，从计数和时间复杂度角度看，第二种方法效率不如第一种方法，但更加准确，能够应对病变严重的眼底图像配准(血管不明显)。但提升的效果需要效率权衡，具体视实验需要效率还是精准度而定。

s108：以获取的眼底荧光造影图像a2作为参考图像，根据仿射变换参数将眼底彩色图像a1和眼底荧光造影图像a2进行配准，具体实施步骤如下：

(1)输入眼底彩色图像a1，将其某个像素点代入fang_final，若变换后的像素点超过眼底荧光造影图像的范围，则抛弃该像素点，若在范围内则保存，最后得到一组像素点，将这组像素点映射到眼底荧光造影图像的像素点数组中，即可完成配准：

arr1_i＝fang_final*arr1_a

其中x1＝x2且y1＝y2。

其中，arr1_a、arr2_a来自s100，可将新的arr_new数组保存为图像n，配准后的图像如图6所示。

s109：对配准后的图像n进行视盘标注，并以视盘中心为圆心，视盘直径dd为基准，画出直径为2dd和3dd的圆，2dd到3dd中的圆构成圆环，如图7所示；若需要测量某一根血管的血管壁厚度，则需要用同样的颜色在圆环内标注该血管中的某个点，存储标注后的图像为q；若需要测量血管的平均血管壁厚度，则无需标注点，存储图像为q0；

s110：根据图像q或q0分别求出图像q中某一根血管的血管壁厚度或图像q0中血管的平均厚度，具体情况如下：

求图像q中某一根血管的血管壁厚度：

(1)输入图像q，读取分辨率即颜色信息，包括宽度x，高度y，以及每个像素点的rgb值，rgb分量值落在0-255区间，存储这些值到数组arr中。由图像q中的标注颜色，得到标注圆环的rgb值b，如红色标注的rgb值为[255,0,0]，绿色标注的rgb值为[0,255,0]。遍历arr数组，读取arr[x][y]的rgb值c，若c等于b，则不处理，若c不等于b则置为[0,0,0]，此值代表黑色：

c＝[r,g,b]

并将此数组保存在内存中。

(2)读取眼底荧光造影图像即眼底造影图像的血管分割图像f1，同样保存其长宽值即像素点颜色分量到数组中，遍历此数组的每个像素点，找到标注为需要测量血管的像素点，以该点为中心，向邻居像素点扩散，每扩散一个点，计数值count1增加一，直到扩散的所有像素点均在圆环之外时结束。此时计数值count1的值则为该眼底荧光造影图像中该血管的像素点总数。用同样的方法处理输入的经过仿射变换后的眼底彩色图像的血管分割图像f2，得到血管像素点总值count2。

(3)用count1减去count2的值即为差值dif，由于dif包含了血管的两侧，所以dif除以2则为血管一侧的差值，用血管一侧的差值除以血管的长度，此长度约等于0.5个dd，得到的结果即为该血管的血管壁厚度估值p：

图像q0中血管的平均厚度：

数组中的点在圆环中时，计数就增加一。这样分别得到的血管内像素点计数为count_total1及count_total2，由于总数包含了8根血管的像素点(四根动脉、四根静脉)。因此平均血管壁厚度为：

尽管上面对本发明说明性的具体实施方式进行了描述，以便于本技术领域的技术人员理解本发明，但应该清楚，本发明不限于具体实施方式的范围，对本技术领域的普通技术人员来讲，只要各种变化在所附的权利要求限定和确定的本发明的精神和范围内，这些变化是显而易见的，一切利用本发明构思的发明创造均在保护之列。