一种基于频域OCT的微血管造影方法与流程

本发明涉及光学检测和生物医学工程领域，特别是涉及一种基于频域光学相关层析(oct)系统进行功能成像的微血管造影方法。

【背景技术】

光学相干层析技术(opticalcoherencetomography，oct)是近年来迅速发展起来的一种成像技术。它在保证一定探测深度的同时，达到了高于其他方法的分辨率和灵敏度，从而以其高分辨率、高灵敏度、非侵入、非接触的特点，在诸多成像方法中占据了一席之地，并有着很好的发展前景，也受到了人们越来越多的关注。光学相干层析成像技术建立在光学、电子学、计算机技术学科的基础上，将半导体和超快激光技术、超灵敏探测、精密自动控制和计算机图像处理等多项技术结合为一整体，是继电子计算机x射线断层扫描成像、超声波成像和核磁共振成像技术之后，又一新的断层成像技术。oct有两大发展趋势，其一是从时域oct到频域oct。时域oct依靠机械方法使参考臂反射镜前后移动来实现深度方向的扫描，而频域oct使用光谱仪分光加傅里叶分析的方法实现深度方向的扫描，成像速度大大提高。其二是从结构成像到功能成像，例如偏振oct、多普勒oct，利用偏振、相位提供更多维度的生物组织信息。

理想的微脉管(直径小于0.3mm的微血管)血流成像技术必须满足几个条件：可以探测到组织内部特定深度的微血液循环；成像响应速率高于血流速度的变化率；能够快速检测血流速度的变化；技术本身应该是非侵入、方便、有效、造价低廉的活体成像术。现有的血流成像方法如彩色超声多普勒技术、荧光血管造影术、激光多普勒测速技术等都不能全部满足这些苛刻要求，它们无法提供较高的空间分辨率，往往还需要静脉注射造影剂来辅助成像。相位分辨多普勒oct技术可以满足这些要求，但缺点是对环境的相位稳定性要求较高，不适合临床应用。

以上背景技术内容的公开仅用于辅助理解本发明的发明构思及技术方案，其并不必然属于本专利申请的现有技术，在没有明确的证据表明上述内容在本专利申请的申请日已经公开的情况下，上述背景技术不应当用于评价本申请的新颖性和创造性。

技术实现要素：

本发明所要解决的技术问题是：弥补上述现有技术的不足，提出一种基于频域oct的微血管造影方法，能检测到更丰富的微血管信息，对检测环境的要求较低，系统组装和调节也较为简化。

本发明的技术问题通过以下的技术方案予以解决：

一种基于频域oct的微血管造影方法，包括以下步骤：s1，构建频域oct系统，所述频域oct系统包括光源、光纤耦合器、参考臂组件、样品臂组件、光谱仪和线阵ccd；所述样品臂组件包括x方向的激光扫描振镜和y方向的激光扫描振镜；所述光源为sled宽带光源；所述样品臂组件的准直光照射在x方向的振镜的转动中心，且偏离y方向的振镜的转动中心2～3mm；通过所述频域oct系统对待测样品进行三维扫描，得到一组x扫描方向的二维光谱；s2，对各张x扫描方向的二维光谱取同一位置的一条线，重构得到一张y扫描方向的二维光谱；将各个位置的线均进行重构，从而得到一组y扫描方向的二维光谱；s3，对各张y扫描方向的二维光谱进行如下处理，得到各y扫描方向的二维光谱对应的血流信息：s31，将当前y扫描方向的二维光谱看成一个二维数组，对该二维数组中的每一列进行插值、减直流、色散补偿的oct数据处理；s32，将色散补偿之后的二维数组在行方向上进行希尔伯特变换构造二维数组的解析信号，对所述二维数组的解析信号在列方向进行fft处理得到新的二维数组，获得待测样品的流速平面强度图和静态平面强度图；s33，利用所述流速平面强度图与静态平面强度图的差值去除所述流速平面强度图中的噪声，得到去噪后的流速平面强度图；然后对所述去噪后的流速平面强度图加伪彩色处理并进行滤波处理，得到当前y扫描方向的二维光谱对应的血流信息；s4，将步骤s3得到的各y扫描方向的二维光谱对应的血管信息拼接成三维血管信息，利用三维可视化软件将三维血管信息转换为三维血管图像。

本发明与现有技术对比的有益效果是：

本发明的血管造影方法，沿x方向进行扫描，但沿y方向进行频移引入及数据处理，这是考虑到系统流体检测的灵敏度和振镜速度密切相关。振镜速度越慢，灵敏度越高，可以检测到更慢的流体。在频域oct系统的三维扫描中，y振镜的速度(一般小于1mm/s)远小于x振镜(一般在300mm/s到1000mm/s)，所以本发明的扫描方式在没有牺牲成像速度的前提下提高了灵敏度，可以检测到速度为几十微米/秒的血液流动，因此能获取更丰富的微血管信息。此外，由于该方法中引入频移时设置了一个速度检测阈值(“偏离y方向的振镜的转动中心)，样品臂准直光偏离y振镜中心的距离越大，速度阈值越高，因此只感应大于该阈值速度的成分，所以对于待测样品的抖动有着更高的容忍度，对检测环境的要求较低。本发明通过样品臂准直光偏离y振镜转动中心的方式引入频移，相对传统的压电陶瓷方法能够更方便地引入频移，对成像速度、扫描长度等参数没有任何制约，系统的组装和调节也更为简化。本发明的微血管造影方法有着极高的分辨率和灵敏度，可以实时、无损的进行活体生物组织的微血管造影，且对环境的要求较低，为疾病评估和生物特征识别提供了丰富的信息。

【附图说明】

图1是本发明具体实施方式中待测成像区域的示意图；

图2是本发明具体实施方式的微血管造影方法中获得的流速平面图像与静态平面图像的示意图；

图3是本发明具体实施方式中提取的二维微血管图；

图4是本发明具体实施方式中获得的三维微血管网络图；

图5是本发明具体实施方式中的频域oct系统的结构示意图。

【具体实施方式】

下面结合具体实施方式并对照附图对本发明做进一步详细说明。

本发明中将oct功能成像应用于微血管造影中，实现对人体活体微循环进行观察研究，无创伤性观察活体微循环，对于疾病的辅助诊断、病情监测、治疗效果的判断、预防保健等具有重大意义。

本具体实施方式中提供一种基于频域oct的微血管造影方法，包括以下步骤：

(1)搭建一套频域oct系统，频域oct系统包括光源、光纤耦合器、参考臂组件、样品臂组件、光谱仪和线阵ccd；所述样品臂组件包括x方向的激光扫描振镜和y方向的激光扫描振镜；所述光源为sled宽带光源；样品臂组件的准直光打在x方向的扫描振镜的转动中心并尽量偏离y方向(偏离2～3mm)的扫描振镜的转动中心以引入固定的载频。通过所述频域oct系统对样品进行三维扫描，得到一组x扫描方向的二维光谱。

(2)每张二维光谱取对应同一位置的一条线，重构得到一张y扫描方向的二维光谱；将各个位置的线均进行重构，从而得到一组y扫描方向的二维光谱。后续独立处理每一张y扫描方向的二维光谱；

(3)y扫描方向的二维光谱可以看成一个二维数组，对该二维数组中的每一列执行oct基本数据处理流程，即插值、减直流、色散补偿；

(4)将色散补偿之后的二维数组先在扫描方向(行方向)做希尔伯特变换构造解析信号，再在波数方向(列方向)做fft(快速傅里叶变换)处理，同时获得流速平面和静态平面的信息；

(5)进行图像处理提取二维血管信息：流速平面与静止平面对应点相减，再加伪彩色使强度和颜色一一对应以增强对比度，最后应用滤波算法进行图像去噪；

(6)对每一幅y扫描方向的二维光谱重复步骤(3)到(5)，拼接成三维血管信息，再利用三维可视化软件完成三维血管网络的可视化。

其中，所述步骤(1)中的频域oct系统使用sled宽带光源，从而系统具有较高的分辨率。优选地，使用中心波长1310nm、半高宽100nm的sled宽带光源，采用像素数为2048的红外线阵ccd接收光谱仪分光后的光谱。由于oct系统的成像深度和ccd像素数成正比，和光源中心波长的平方成正比。光源中心波长一定的情况下，光源半高宽越高，纵向分辨率越高。因此在可选的范围内取上述像素数较多的ccd，半高宽高达100nm的光源，从而系统的纵向分辨率较高。

所述步骤(1)中的样品臂准直光束直径为2.8mm，x振镜尺寸12.5mm*12.5mm，y振镜尺寸12.5mm*17.5mm。在其他条件一定的情况下，样品臂准直光束直径越小，流体检测的灵敏度越高(即可以检测到更小速度的流动)。由于系统调节时需要使光束偏离y振镜中心，偏移量越大对应的速度阈值越大。所以y振镜的尺寸越大(12.5mm*17.5mm)，阈值速度的选择范围就越大，方便系统的调节。

所述步骤(1)中的样品例如内有散射体流淌的毛细玻璃管、人的手指等等，样品可辅助胶带、机械件进行较严格的固定。系统可感应沿光束方向的速度分量。

所述步骤(1)中的三维扫描由x、y两个激光扫描振镜完成，一般将扫描速度快的方向称为x方向(x方向的振镜的速度范围一般在300mm/s～1000mm/s)，扫描速度慢的方向称为y方向(y方向的振镜的速度一般小于1mm/s)。x振镜快速往复运动做二维扫描，y振镜缓慢步进运动做三维扫描。一组实验参数中，x方向振镜的速度为600mm/s，y方向振镜的速度为0.83mm/s。x方向扫描密度为60lines/mm，y方向扫描密度为600lines/mm。

所述步骤(2)中，“对应”位置即每张二维光谱图的同一位置。例如步骤(1)得到一组共300张x扫描方向的二维光谱图，每张二维光谱图中有500个一维光谱(即500条线)。如果要研究每张二维光谱中第一条线对应的切片，就取每张二维光谱图中的第一条线，300个二维光谱一共可以取到300条线，重构出来的就是一个y扫描方向的二维光谱，该二维光谱中有300条线。这一组300张x扫描方向的二维光谱图，最终可重构得到一组共500张y扫描方向的二维光谱(每张光谱中有300条线)。

所述步骤(3)中的插值是指在调试系统时通过反射镜实验得到波数k随ccd像元均匀变化的关系(k～n)，然后利用这一关系将光强随ccd像元均匀变化的函数(i～n)转换为光强随波数k均匀变化的函数(i～k)。此处的波数k和n的关系是频域oct系统中光谱仪这一硬件决定的，是光谱仪的硬件性质，这个关系反映ccd的不同位置对应的波数是多少。频域oct最大的特色是可以获取样品的深度断层信息，所以oct数据处理的核心步骤是对(i～k)关系做fft以解析出光强随深度变化的关系。然而光谱仪直接得到的光谱是(i～n)关系，所以通过该插值的步骤将(i～n)关系转化为(i～k)关系。

所述步骤(3)中的减直流是为了减去干涉光谱中的直流量以减小零光程附近的噪声。方法是将该二维图的光谱平均值做为该二维图的直流量从当前y扫描方向的二维光谱中的每条线中减去。例如，一张y扫描方向的二维数组有300列，通常称这个二维数组有300条线。每一条线都是单独进行处理的。这里的“二维光谱的光谱平均值”相当于是这300条线相加除以300作为直流量，然后从每条线都减去这个直流量。

所述步骤(3)中的色散补偿是指在当前y扫描方向的二维光谱的每条线的相位项中用软件减去色散附加相位，即关于k的二阶及以上阶数的高阶项。

此处步骤(3)中的插值、减直流、色散补偿等处理，均是进行fft处理之前的预处理操作。因为获取的光谱中叠加了色散附加相位，如果直接做fft会引入噪声，影响图像的纵向分辨率。所以通过色散补偿，减去色散附加相位。

所述步骤(4)中，以处理大小为2048*1200的二维数组为例，即红外线阵ccd有2048个像素，二维数组有1200条线。首先单独处理该二维数组的每一行，即单独处理长度为1200的一维数组。对一维数组进行希尔伯特变换；以原一维数组为实部，变换后的一维数组为虚部，构造一个新的复数一维数组，定义该复数一维数组为原一维数组的解析信号。将以上操作重复2048次，就得到了一个新的二维数组，该二维数组的每一行都是原二维数组对应行的解析信号。

然后单独处理这个新的二维数组的每一列，即单独处理长度为2048的一维数组。对一维数组进行fft，重复1200次，就得到了一个新的复数二维数组。对这个新的复数二维数组取模，从而同时得到流速平面强度图和静态平面强度图。需说明的是，原始数据，即fft之前的二维数组的y轴是波数k，数组中数字的物理意义是光谱值。进行fft处理后的二维数组的y轴是深度z，数组中数字的物理意义是强度值。光谱值转换为强度值，即由k空间转换到了z空间(深度空间)。步骤s32后则是对强度值进行处理。

上述过程中，进行希尔伯特变换构造的解析信号是一个复数信号，因此解析信号的fft是一个非对称的双边谱，样品流动部分和静止部分相位符号的差异导致它们分别成像在不同的半平面上。流速平面强度图和静态平面强度图分别为经过fft处理后得到的强度图的正负平面。

具体地，fft处理后得到的新的复数二维数组取模后，得到了一个大小为2048*1200的实数二维数组。这个二维数组可以横向切开，分成两个1024*1200的二维数组。其中，上方的二维数组为fft得到的正频率成分，下方的二维数组为fft得到的负频率成分。这两个二维数组，一个对应流速平面强度图，另一个对应静态平面强度图。具体是上方的二维数组对应流速平面还是下方的二维数组对应流速平面取决于样品臂光束偏移y振镜中心的方向，对最终的造影结果没有影响。在微血管造影的应用中，所有流速部位的光强值之和是小于所有静态部位的光强值之和。因此，分别计算两个二维数组的光强值之和，和较大的二维数组是静态平面强度图，和较小的二维数组流速平面强度图。

步骤(4)中fft变换后得到的新的复数二维数组后，可保留所有的数据，也可适当截取部分数据，都不影响后续的操作和微血管造影效果。由于oct的成像深度和穿透深度有限，fft后的高频成分没有参考意义。因此，优选地，可保留正频率成分的前500个点和负频率成分的前500个点，共计1000个点进行后续处理。即将原来的2048*1200的复数二维数组截取出最前面的部分行复数(例如最前面的500行复数)，以及最后面的部分行复数(同样是500行复数)，由这个1000*1200的复数二维数组进行后续的取模、横向切开后分别获得对应的流速平面强度图和静态平面强度图。

所述步骤(5)中相减指流速平面强度图对应点减去静态平面强度图对应点。利用静态平面图像，对流速平面图像做相减处理，得到新的流速平面图像。若小于零，说明是没消干净的镜像，将流速平面的该点的强度值置为所述流速平面强度图中各点对应的强度值中的最小强度值；若大于零，说明是流速像，保持该点的原强度值。本具体实施方式中，目的是获得流速信息；静态平面只是一个参考面，用来做相减处理以判断是否应该将流速平面的对应点置为黑色。如果相减后的强度小于零，下一个步骤就将该点(流速平面)的颜色设置为黑色(取各点对应的强度值中的最小强度值)。如果大于零，就保留流速平面强度图中这一点的原强度值。

步骤(5)相减处理后，加伪彩色使强度和颜色一一对应以增强对比度时，“加伪彩色”是将前述fft处理后的强度值转换为颜色值(rgb)。综上，变换的步骤依次是光谱值-强度值-颜色值(rgb)。

所述步骤(5)中运用bm3d滤波算法对新的流速平面图像进行滤波去噪处理，目的是从新的流速平面图像当中减弱加性的高斯白噪声。bm3d噪声的标准差取为50。

所述步骤(6)中的可视化需先将一系列二维图转化为一个.raw文件，再用inviwo调节传递函数，赋予不同灰度以不同的颜色。

具体地，本具体实施方式中，搭建的频域oct系统的结构示意图如图5所示，样品臂准直光束直径为2.8mm，x振镜尺寸12.5mm*12.5mm，y振镜尺寸12.5mm*17.5mm。使样品臂准直光打在x振镜的转动中心并尽量偏离y振镜的转动中心以引入固定的载频。将手指放在样品台上，对图1条形框所示的区域进行三维扫描，得到一组x扫描方向的二维光谱。该区域大小为1mm*2mm，x方向扫描密度为60lines/mm，y方向扫描密度为600lines/mm。

每张二维光谱取对应位置的一条线，重构为60幅y扫描方向的二维光谱。后续独立处理每一幅y扫描方向的二维光谱，这个二维光谱可以看成一个尺寸为1200*2048的二维数组，对该二维数组中的每一列执行oct基本数据处理流程，即插值、减直流、色散补偿。将色散补偿之后的二维数组先在扫描方向(行方向)做希尔伯特变换构造解析信号，再在波数方向(列方向)做fft，同时获得流速平面和静态平面的信息如图2所示，箭头指向的部分为微血管。

进行图像处理提取二维血管信息。流速平面对应点减去静态平面对应点。若小于零，说明是没消干净的镜像，将该点置黑；若大于零，说明是流速像，保持该点的原值。然后应用bm3d滤波算法进行图像去噪，去噪后的结果如图3所示。

对每一幅y扫描方向的二维光谱重复以上步骤，拼接成三维血管信息，再利用三维可视化软件inviwo调节传递函数，赋予不同灰度以不同的颜色，完成三维血管网络的可视化如图4所示，探测到的微血管直径约为80um。

本具体实施方式的微血管造影方法，相对于其他的微血管造影方法如fa(荧光素血管造影)、icga(吲哚青绿血管造影)，基于oct的微血管造影方法不需要注射造影剂，是非侵入的检测方法，更加安全快捷；oct方法本身也具有十几微米的横向分辨率和纵向分辨率，因此可以获得更高分辨率的血管图像。此外，oct具有深度分辨能力，因此可以查看任意深度的血管切片图，更容易掌握微血管的三维结构信息。

相对其他的oct微血管造影方法，本具体实施方式的微血管造影方法也有如下优势。本方法通过样品臂准直光偏离y振镜转动中心的方式引入频移，相对传统的压电陶瓷方法能够更方便地引入频移，对成像速度、扫描长度等参数没有任何制约，系统的组装和调节也更为简化。

和相位分辨多普勒oct相比，本具体实施方式的微血管造影方法通过引入频移设置了一个速度检测阈值，样品臂准直光偏离y振镜中心的距离越大，速度阈值越高，因此只感应大于该阈值速度的成分，所以对于待测样品的抖动有着更高的容忍度，对检测环境的要求较低。

在扫描方式上，本具体实施方式的微血管造影方法沿x方向进行扫描，但沿y方向进行频移引入及数据处理，这是由于系统流体检测的灵敏度和振镜速度密切相关。振镜速度越慢，灵敏度越高，可以检测到更慢的流体。在三维扫描中，y振镜的速度远小于x振镜，所以该扫描方式在没有牺牲成像速度的前提下提高了灵敏度，可以检测到几十微米的血液流动，因此能获取更丰富的微血管信息。

本具体实施方式的微血管造影方法有着极高的分辨率和灵敏度，可以实时、无损的进行活体生物组织的微血管造影，且对环境的要求较低，为疾病评估和生物特征识别提供了丰富的信息。

以上内容是结合具体的优选实施方式对本发明所作的进一步详细说明，不能认定本发明的具体实施只局限于这些说明。对于本发明所属技术领域的普通技术人员来说，在不脱离本发明构思的前提下做出若干替代或明显变型，而且性能或用途相同，都应当视为属于本发明的保护范围。