CN114894794B - 干涉图的消除伪影的方法 - Google Patents

Abstract

本公开描述了一种干涉图的消除伪影的方法。该方法应用于利用平面干涉仪对鲜活样本的成像中，包括：获取多张干涉图并组成第一干涉图像集，多张干涉图是通过平面干涉仪针对鲜活样本进行采集而获得；对第一干涉图像集进行奇异值分解以获取时间特征矩阵和奇异值矩阵，时间特征矩阵的每列为时间特征向量；计算时间特征向量的波动程度；基于时间特征向量的波动程度更新奇异值矩阵中的奇异值；并且基于更新后的奇异值矩阵重建第一干涉图像集并作为第二干涉图像集以基于第二干涉图像集生成彩色图像，彩色图像反映鲜活样本的内部动态信号。由此，能够通过非侵入式的方式获得生物组织的内部动态信号。

Description

干涉图的消除伪影的方法

本申请是申请日为2021年04月08日、申请号为2021103780927、发明名称为针对平面干涉仪的信号处理的成像系统及成像方法的分案申请。

技术领域

本公开大体涉及一种干涉图的消除伪影的方法。

背景技术

生物组织的内部动态信号能够为生物医学领域的科学研究及临床诊疗提供重要参考，例如在组织病理学或干细胞组织工程学领域，生物组织的内部动态信号的强弱常常用于细胞类型区分或细胞培养检测。因此，如何获取生物组织的内部动态信号已经成为生物医学领域热门的研究方向。目前可以通过侵入式的方法(例如基于荧光染色的方法)获取生物组织的内部动态信号。

然而，侵入式的方法一般具有破坏性，不适用于鲜活生物组织的观察和测量，而通过非侵入式的方法(例如利用基于双光束等厚干涉原理的平面干涉仪对生物组织进行测量)获得的生物组织的结构成像却存在无法反映动态信号的问题，无法实现生物组织的内部动态信号的测量。

发明内容

本公开是鉴于上述的状况而提出的，其目的在于提供一种能够通过非侵入式的方式获得生物组织的内部动态信号的针对平面干涉仪的信号处理的成像系统及成像方法。

为此，本公开的第一方面提供了一种针对平面干涉仪的信号处理的成像系统，所述平面干涉仪基于双光束等厚干涉原理进行测量，所述成像系统包括：获取模块、消除伪影模块、信噪比增强模块和动态成像模块；所述获取模块用于获取利用所述平面干涉仪对鲜活样本进行连续采集获得的多张干涉图，所述多张干涉图按照采集时间组成第一干涉图像集；所述消除伪影模块用于对所述第一干涉图像集进行奇异值分解以获取时间特征矩阵和奇异值矩阵，所述时间特征矩阵的每列为时间特征向量，计算所述时间特征向量的波动程度并基于所述波动程度筛选符合阈值条件的时间特征向量以将该时间特征向量对应的所述奇异值矩阵中的奇异值更新为预设奇异值，然后基于更新后的奇异值矩阵重建所述第一干涉图像集并作为第二干涉图像集；所述信噪比增强模块利用预设长度的滑动窗对所述第二干涉图像集中的各个滑动窗中的像素值进行信噪比增强处理以获取第三干涉图像集，其中，所述滑动窗按照预设步长沿着所述采集时间的方向移动；以及所述动态成像模块用于获取所述第三干涉图像集中的多张干涉图的相同位置的多个像素值以形成按照所述采集时间排序的像素序列，基于各个位置的像素序列的时域信息和频域信息获取各个位置的三维颜色空间的分量并基于所述三维颜色空间的分量将所述第三干涉图像集转为彩色图像，其中，所述彩色图像反映所述鲜活样本的内部动态信号，所述三维颜色空间的分量包括色调分量、饱和度分量和亮度分量，基于所述频域信息获取目标频率，所述色调分量与所述目标频率对应。在本公开中，成像系统基于奇异值分解获得时间特征向量以通过时间特征向量的波动程度对多张干涉图的伪影进行消除并利用滑动窗对多张干涉图进行信噪比增强处理，最终基于多张干涉图的时域信息和频域信息生成彩色图像。由此，能够降低伪影或者噪声对彩色图像的图像质量的影响，有利于基于彩色图像进行进一步地研究，另外通过本公开生成的彩色图像能够反映生物组织的内部动态信号并且获得该彩色图像的方式是非侵入式。

另外，在本公开的第一方面所涉及的成像系统中，可选地，所述成像系统还包括标定模块，所述标定模块用于建立所述目标频率与颜色的线性对应关系以对所述目标频率与所述颜色的对应关系进行标定。由此，能够直观地获得彩色图像中的颜色对应的频率。

另外，在本公开的第一方面所涉及的成像系统中，可选地，所述消除伪影模块利用所述时间特征向量的累积过零率表示所述时间特征向量的波动程度，第i列时间特征向量的累积过零率表示为：D_ZRC_i＝|ZRC_i+1-ZRC_i|，其中，i为所述时间特征矩阵的列索引，ZRC_i为第i列时间特征向量的过零率。由此，能够利用累积过零表示时间特征向量的波动程度。

另外，在本公开的第一方面所涉及的成像系统中，可选地，所述阈值条件为各个时间特征向量的累积过零率大于预设值，其中，所述预设值为3倍的所述累积过零率的标准差。由此，能够通过消除累积过零率的离散程度达到预设值程度的信号以消除伪影。

另外，在本公开的第一方面所涉及的成像系统中，可选地，所述预设长度大于1，所述滑动窗的维度与所述第二干涉图像集的维度一致。

另外，在本公开的第一方面所涉及的成像系统中，可选地，在所述信噪比增强处理中，获取各个滑动窗中的像素值的平均值，对各个滑动窗中的像素值与所述平均值之间的差值求累积和以获取累积值，将所述累积值的绝对值与所述预设长度相除以获取平均累积值，将所述平均累积值作为所述第三干涉图像集的像素值。由此，能够基于滑动窗和累积和的算法对第二干涉图像集进行信噪比增强处理。

另外，在本公开的第一方面所涉及的成像系统中，可选地，所述动态成像模块基于各个位置的像素序列获取功率谱密度并将所述功率谱密度的带宽的倒数作为所述饱和度分量；所述动态成像模块对各个位置的像素序列进行傅立叶变换以获取频率序列，基于所述频率序列和所述功率谱密度获取所述目标频率并作为所述色调分量；所述动态成像模块基于各个位置的像素序列获取该像素序列的标准差或方差并将该标准差或该方差作为所述亮度分量。由此，能够基于各个位置的像素序列获得饱和度分量、色调分量和亮度分量。

本公开的第二方面提供了一种针对平面干涉仪的信号处理的成像方法，所述平面干涉仪基于双光束等厚干涉原理进行测量，所述成像方法包括：获取利用所述平面干涉仪对鲜活样本进行连续采集获得的多张干涉图，所述多张干涉图按照采集时间组成第一干涉图像集；对所述第一干涉图像集进行奇异值分解以获取时间特征矩阵和奇异值矩阵，所述时间特征矩阵的每列为时间特征向量，计算所述时间特征向量的波动程度并基于所述波动程度筛选符合阈值条件的时间特征向量以将该时间特征向量对应的所述奇异值矩阵中的奇异值更新为预设奇异值，然后基于更新后的奇异值矩阵重建所述第一干涉图像集并作为第二干涉图像集；利用预设长度的滑动窗对所述第二干涉图像集中的各个滑动窗中的像素值进行信噪比增强处理以获取第三干涉图像集，其中，所述滑动窗按照预设步长沿着所述采集时间的方向移动；并且获取所述第三干涉图像集中的多张干涉图的相同位置的多个像素值以形成按照所述采集时间排序的像素序列，基于各个位置的像素序列的时域信息和频域信息获取各个位置的三维颜色空间的分量并基于所述三维颜色空间的分量将所述第三干涉图像集转为彩色图像，其中，所述彩色图像反映所述鲜活样本的内部动态信号，所述三维颜色空间的分量包括色调分量、饱和度分量和亮度分量，基于所述频域信息获取目标频率，所述色调分量与所述目标频率对应。在本公开中，成像方法基于奇异值分解获得时间特征向量以通过时间特征向量的波动程度对多张干涉图的伪影进行消除并利用滑动窗对多张干涉图进行信噪比增强处理，最终基于多张干涉图的时域信息和频域信息生成彩色图像。由此，能够降低伪影或者噪声对彩色图像的图像质量的影响，有利于基于彩色图像进行进一步地研究，另外通过本公开生成的彩色图像能够反映生物组织的内部动态信号并且获得该彩色图像的方式是非侵入式。

另外，在本公开的第二方面所涉及的成像方法中，可选地，建立所述目标频率与颜色的线性对应关系以对所述目标频率与所述颜色的对应关系进行标定。由此，能够直观地获得彩色图像中的颜色对应的频率。

另外，在本公开的第二方面所涉及的成像方法中，可选地，利用所述时间特征向量的累积过零率表示所述时间特征向量的波动程度，第i列时间特征向量的累积过零率表示为：D_ZRC_i＝|ZRC_i+1-ZRC_i|，其中，i为所述时间特征矩阵的列索引，ZRC_i为第i列时间特征向量的过零率；所述阈值条件为各个时间特征向量的累积过零率大于预设值，其中，所述预设值为3倍的所述累积过零率的标准差。由此，能够通过消除累积过零率的离散程度达到预设值程度的信号以消除伪影。

根据本公开，能够提供一种通过非侵入式的方式获得生物组织的内部动态信号的针对平面干涉仪的信号处理的成像系统及成像方法。

附图说明

现在将仅通过参考附图的例子进一步详细地解释本公开，其中：

图1是示出了本公开示例所涉及的针对平面干涉仪的信号处理的成像方法的应用场景图示意图。

图2是示出了本公开示例所涉及的针对平面干涉仪的信号处理的成像系统的框图。

图3是示出了本公开示例所涉及的干涉图的示意图。

图4(a)是示出了本公开示例所涉及的为未消除伪影情况下生成的彩色图像的示意图。

图4(b)是示出了本公开示例所涉及的为消除伪影情况下生成的彩色图像的示意图。

图4(c)是示出了本公开示例所涉及的消除的伪影的示意图。

图5是示出了本公开示例所涉及的滑动窗移动的示意图。

图6是示出了本公开示例所涉及的带有标定信息的彩色图像的示意图。

图7是示出了本公开示例所涉及的针对平面干涉仪的信号处理的成像方法的流程图。

具体实施方式

以下，参考附图，详细地说明本公开的优选实施方式。在下面的说明中，对于相同的部件赋予相同的符号，省略重复的说明。另外，附图只是示意性的图，部件相互之间的尺寸的比例或者部件的形状等可以与实际的不同。需要说明的是，本公开中的术语“包括”和“具有”以及它们的任何变形，例如所包括或所具有的一系列步骤或单元的过程、方法、系统、产品或设备不必限于清楚地列出的那些步骤或单元，而是可以包括或具有没有清楚地列出的或对于这些过程、方法、产品或设备固有的其它步骤或单元。本公开所描述的所有方法可以以任何合适的顺序执行，除非在此另有指示或者与上下文明显矛盾。

本公开涉及的针对平面干涉仪的信号处理的成像系统及成像方法能够基于由平面干涉仪采集的多张干涉图通过非侵入式的方式获得反映生物组织的内部动态信号的彩色图像。本公开涉及的成像方法应用于成像系统(稍后描述)。以下结合附图进行详细描述本公开。此外，本公开的示例描述的应用场景是为了更加清楚的说明本公开的技术方案，并不构成对于本公开提供的技术方案的限定。

图1是示出了本公开示例所涉及的针对平面干涉仪的信号处理的成像方法的应用场景图示意图。如图1所示，可以利用平面干涉仪110对生物组织120进行测量以获得多张干涉图130。在完成多张干涉图130采集后，可以将多张干涉图130提交至服务器140，服务器140可以通过执行存储在服务器140上的计算机程序指令以实现成像方法，通过该成像方法接收多张干涉图130并生成彩色图像150。

在一些示例中，平面干涉仪110可以基于双光束等厚干涉原理进行测量。在一些示例中，平面干涉仪110可以是二维平面干涉仪，二维平面干涉仪可以是由2D相机进行干涉图130的采集。在一些示例中，获得多张干涉图130的设备可以是任意基于双光束等厚干涉原理进行测量的设备。在一些示例中，生物组织120可以是介于细胞及器官之间的细胞架构。

在一些示例中，服务器140可以包括一个或多个处理器和一个或多个存储器。其中，处理器可以包括中央处理单元、图形处理单元以及能够处理数据的其它任何电子部件，能够执行计算机程序指令。存储器可以用于存储计算机程序指令。在一些示例中，服务器140通过执行存储器上的计算机程序指令可以实现成像方法。在一些示例中，服务器140也可以是云端服务器。

以下，结合附图详细描述本公开涉及的成像系统1。本公开涉及的成像系统1用于实施上述的成像方法。图2是示出了本公开示例所涉及的针对平面干涉仪的信号处理的成像系统1的框图。图3是示出了本公开示例所涉及的干涉图的示意图。

如图2所示，在一些示例中，成像系统1可以包括获取模块10，获取模块10可以用于获取多张干涉图并组成第一干涉图像集。干涉图可以为波的干涉所形成的图样。在一些示例中，多张干涉图的尺寸可以一致。例如多张干涉图的尺寸可以均为1440×1440或1024×1024。但本公开的示例不限于此，在另一些示例中，干涉图的尺寸可以由平面干涉仪中用于采集干涉图的相机所决定。在一些示例中，多张干涉图的数量可以大于等于预设数量。例如预设数量可以为16，多张干涉图的数量例如可以为16、32、48、64、100、512或1024等。作为干涉图的示例，图3示出了由平面干涉仪测量获得的干涉图的示意图。

另外，在一些示例中，可以利用平面干涉仪对样本进行采集以获得多张干涉图。在一些示例中，样本可以是鲜活样本。例如小鼠视网膜的鲜活样本。在一些示例中，鲜活样本可以包括内容物。在一些示例中，内容物可以为细胞。在一些示例中，多张干涉图可以是通过连续采集获得的。在一些示例中，可以基于选择的曝光频率(也即每秒曝光的次数)进行连续采集。例如可以选择曝光频率为100次/秒连续采集以获得512张干涉图，在这种情况下，完成512张干涉图的采集需耗时5.12秒。

另外，在一些示例中，可以将多张干涉图按照采集时间组成第一干涉图像集。由此，能够方便后续同时对多张干涉图进行处理。

图4(a)是示出了本公开示例所涉及的为未消除伪影情况下生成的彩色图像的示意图。图4(b)是示出了本公开示例所涉及的为消除伪影情况下生成的彩色图像的示意图。图4(c)是示出了本公开示例所涉及的消除的伪影的示意图。

如图2所示，在一些示例中，成像系统1可以包括消除伪影模块11，消除伪影模块11可以用于对第一干涉图像集进行消除伪影处理以获取第二干涉图像集。一般而言，伪影(Artifacts)是指原本被扫描物体例如鲜活样本并不存在而在最终成像的图像上却出现的各种形态的影像。一般而言，造成伪影原因往往有多种。例如鲜活样本在活跃状态会有运动可造成伪影，这种运动在空间层面非常微小，例如是微米级别的。又例如环境空气的扰动也可能带来振动，进而产生伪影。在这种情况下，对第一干涉图像集中的伪影进行消除能够提高后续生成的彩色图像的图像质量。

作为未消除伪影和消除伪影生成的彩色图像的对比示例。图4(a)和图4(b)分别了示出了小鼠视网膜的鲜活样本的未消除伪影情况下生成的彩色图像和消除伪影情况下生成的彩色图像，图4(c)示出了小鼠视网膜的鲜活样本中消除的伪影。从图4(a)、图4(b)和图4(c)中可以看出，消除了区域B所示的伪影的区域C相对未消除伪影的区域A对比度明显增强。由此，能够识别被伪影遮挡的小细胞。这里区域A、区域B和区域C对应生成的彩色图像中的同一个位置。

在一些示例中，消除伪影模块11可以基于奇异值分解对第一干涉图像集进行消除伪影处理以获取第二干涉图像集。在一些示例中，在消除伪影处理中，可以对第一干涉图像集进行奇异值分解(Singular Value Decomposition，SVD)以获取时间特征矩阵和奇异值矩阵。

具体地，可以将第一干涉图像集转成二维数据S，该二维数据S的大小可以为m×n，其中，m可以为一张干涉图的分辨率(尺寸)，n可以为多张干涉图的数量。在一些示例中，二维数据S中每行可以包括各张干涉图的所有像素值。例如，连续采集512张干涉图，每张干涉图的尺寸为1440×1440，则m可以为1440×1440，n可以为512。

在一些示例中，对二维数据S进行奇异值分解以将二维数据S表示为：

S＝U∑V，

其中，U为空间特征矩阵，∑为奇异值矩阵，V为时间特征矩阵。在一些示例中，时间特征矩阵的大小可以为n×n。例如假设连续采集512张干涉图，则时间特征矩阵的大小可以512×512。奇异值矩阵可以为对角矩阵，对角线上的元素是一个奇异值。在一些示例中，奇异值矩阵的大小可以取n×n。例如假设连续采集512张干涉图，则奇异值矩阵的大小可以512×512。

在一些示例中，在获取时间特征矩阵和奇异值矩阵后，可以将时间特征矩阵的每列作为时间特征向量。由于伪影对应的时间特征向量相对正常的干涉信号对应的时间特征向量有较大的区别，因此，基于时间特征向量的波动程度能够区分出伪影。在一些示例中，可以基于时间特征向量的波动程度更新奇异值矩阵中的奇异值。具体地，可以计算时间特征向量的波动程度，基于波动程度筛选符合阈值条件的时间特征向量，然后将该时间特征向量对应的奇异值矩阵中的奇异值更新为预设奇异值。由此，能够基于时间特征向量的波动程度对奇异值矩阵中的奇异值进行更新。在一些示例中，预设奇异值可以为0。

在一些示例中，波动程度可以利用时间特征向量的累积过零率进行表示。也即，消除伪影模块11可以利用时间特征向量的累积过零率表示时间特征向量的波动程度。如上所述，可以将时间特征矩阵的每列作为时间特征向量。在一些示例中，第i列时间特征向量的累积过零率可以表示为：

D_ZRC_i＝|ZRC_i+1-ZRC_i|，

其中，i为时间特征矩阵的列索引，ZRC_i为第i列时间特征向量的过零率。过零率可以指一个信号的符号变化的比率，在本公开中可以指时间特征向量中的多个值形成的波动图(也可以称为变化曲线图)经过零的次数。由此，能够利用累积过零表示时间特征向量的波动程度。在一些示例中，可以绘制出每个时间特征向量的曲线并基于曲线求过零率。在一些示例中，i可以取值1至n-1，n可以为多张干涉图的数量。

但本公开的示例不限于此，在另一些示例中，波动程度可以利用其他方式进行表示。

在一些示例中，阈值条件可以为各个时间特征向量的累积过零率大于预设值。在一些示例中，预设值可以为3倍的累积过零率的标准差。也即，阈值条件可以表示为：

D_ZRC_i>3×std(D_ZRC)，

其中，D_{_}ZRC_i为第i列时间特征向量的累积过零率，D_{_}ZRC为累积过零率的集合。由此，能够通过消除累积过零率的离散程度达到预设值程度的信号以消除伪影。但是本公开的示例不限于此，在另一些示例中，可以根据实验经验或实验效果设置预设值。

如上所述，可以基于时间特征向量的波动程度更新奇异值矩阵中的奇异值。在一些示例中，基于更新后的奇异值矩阵可以重建第一干涉图像集并作为第二干涉图像集。由此，能够基于时间特征向量的波动程度对第一干涉图像集进行消除伪影处理。具体地，可以将上述的空间特征矩阵、更新后的奇异值矩阵以及时间特征矩阵相乘以获取更新后的二维数据，然后基于更新后的二维数据获取与第一干涉图像集对应的第二干涉图像集。

但是本公开的示例不限于此，在另一些示例中，成像系统1也可以不包括消除伪影模块11。在这种情况下，后续模块例如信噪比增强模块12或动态成像模块13可以对第一干涉图像集进行处理。

图5是示出了本公开示例所涉及的滑动窗移动的示意图。如图2所示，在一些示例中，成像系统1可以包括信噪比增强模块12，信噪比增强模块12可以用于对第二干涉图像集进行信噪比增强处理以获取第三干涉图像集。在一些示例中，信噪比增强模块12可以利用预设长度的滑动窗对第二干涉图像集中的各个滑动窗中的像素值进行信噪比增强处理以获取第三干涉图像集。在一些示例中，预设长度可以大于1。在一些示例中，预设长度可以为8。

在一些示例中，第二干涉图像集可以表示为三维的矩阵，其中，第一维和第二维可以对应干涉图的尺寸，第三维对应采集顺序(也可以称为采集时间)。在一些示例中，滑动窗的维度与第二干涉图像集的维度可以一致。在一些示例中，滑动窗的维度的第一维和第二维可以与干涉图的尺寸一致。在这种情况下，能够使滑动窗基于干涉图的尺寸沿着采集顺序移动。

在一些示例中，滑动窗可以沿着采集时间的方向移动。作为滑动窗移动的示例，图5是示出了预设长度为2的滑动窗H沿着第二干涉图像集L移动的示意图，其中第二干涉图像集L包含多张干涉图，例如多张干涉图的数量可以为n，则多张干涉图可以包括p1、p2、p3、……、pn等。

在一些示例中，滑动窗可以按照预设步长沿着采集时间的方向移动。在一些示例中，预设步长可以为预设长度的预设比例。在一些示例中，预设比例可以小于1。例如预设比例可以为25％、50％、75％、1/3或2/3等。具体地，假设预设比例为50％，第二干涉图像集中的干涉图数量为512，预设长度为8，则预设步长为4，滑动窗的个数为512/4个，也即128个。在这种情况下，滑动窗之间存在重叠，后一个滑动窗能够结合前一个滑动窗中的部分干涉图进行信噪比增强处理。但本公开的示例不限于此，在另一些示例中，可以根据实际情况设置预设步长和预设长度，例如预设长度可以为6，预设步长可以为2。

在一些示例中，在信噪比增强处理中，可以获取各个滑动窗中的像素值的平均值。在一些示例中，可以对各个滑动窗中的像素值与平均值之间的差值求累积和以获取累积值。在一些示例中，可以将累积值的绝对值与预设长度相除以获取平均累积值，将平均累积值作为第三干涉图像集的像素值。由此，能够基于滑动窗和累积和的算法对第二干涉图像集进行信噪比增强处理。

具体地，第三干涉图像集S₃可以表示为：

其中，λ为预设长度(也可以称为滑动窗长度)，W_k为第k个滑动窗内的像素值，为第k个滑动窗中的像素值的平均值，CumSum为累积和函数，K为滑动窗的个数。

但是本公开的示例不限于此，在另一些示例中，成像系统1也可以不包括信噪比增强模块12。在这种情况下，后续模块例如动态成像模块13可以对第一干涉图像集或第二干涉图像集进行处理。

如图2所示，在一些示例中，成像系统1可以包括动态成像模块13，动态成像模块13可以基于第三干涉图像集生成彩色图像。彩色图像可以反映鲜活样本的内部动态信号。在一些示例中，可以基于第三干涉图像集获取三维颜色空间的分量并基于三维颜色空间的分量将第三干涉图像集转为彩色图像。在一些示例中，三维颜色空间的分量可以包括色调分量、饱和度分量和亮度分量(也可以称为明度分量)。在一些示例中，三维颜色空间可以为HSV(Hue:色相，Saturation:饱和度，Value:明度)空间或HSI(Hue:色相，Saturation:饱和度，Intensity:亮度)空间。由此，能够通过HSV空间或HSI空间生成彩色图像。在一些示例中，彩色图像可以为RGB(Red:红，Green:绿，Blue:蓝)空间的图像。

在一些示例中，可以获取第三干涉图像集中的多张干涉图的相同位置的多个像素值以形成按照采集时间排序的像素序列，基于各个位置的像素序列的时域信息和频域信息获取各个位置的三维颜色空间的分量，并基于三维颜色空间的分量将第三干涉图像集转为彩色图像。

一般而言，频域信息可以包括信号的频率分量以及各个频率分量的幅值大小，时域信息可以包括信号的幅值信息，可以测量到时域信号波形随着时间的变化快慢。如上所述，像素序列是各个位置的多个像素值按照采集时间排序的形成的，因此，在本公开中，时域信息可以指像素序列中的随采集时间变化的像素值。在一些示例中，可以通过傅立叶级数或傅立叶变换将信号从时域变换到频域以获取频域信息。

如上所述，三维颜色空间的分量可以包括色调分量、饱和度分量和亮度分量。在一些示例中，动态成像模块13可以基于频域信息获取色调分量。在一些示例中，动态成像模块13可以基于频域信息获取目标频率并基于目标频率获取色调分量。

在一些示例中，色调分量可以与目标频率对应。在一些示例中，可以将目标频率作为色调分量。在一些示例中，目标频率可以为平均频率或中值频率。在一些示例中，动态成像模块13可以对各个位置的像素序列进行傅立叶变换以获取频率序列(也即频域信息)，基于频率序列和功率谱密度(稍后描述)获取平均频率并作为色调分量。在一些示例中，可以对频率序列和功率谱密度求点积以获取平均频率。由此，能够基于各个位置的像素序列获得色调分量。

如上所述，目标频率可以为中值频率。在一些示例中，中值频率可以为将功率谱密度分成上、下两个相等功率面积的频率。在一些示例中，中值频率MDF可以满足公式：

其中，P(f)为功率谱密度，f为频率。

在一些示例中，动态成像模块13可以基于各个位置的像素序列(也即时域信息)获取饱和度分量。在一些示例中，动态成像模块13可以基于各个位置的像素序列获取功率谱密度并基于功率谱密度获取饱和度分量。在一些示例中，可以对功率谱密度进行归一化再基于归一化后的功率谱密度获取饱和度分量。在一些示例中，动态成像模块13可以基于各个位置的像素序列获取功率谱密度并将功率谱密度的带宽的倒数作为饱和度分量。在一些示例中，可以基于各个位置的像素序列并利用Welch法获取功率谱密度。由此，能够基于各个位置的像素序列获得饱和度分量。

在一些示例中，动态成像模块13可以基于各个位置的像素序列(也即时域信息)获取该像素序列的标准差或方差并将该标准差或该方差作为亮度分量。由此，能够基于各个位置的像素序列获得亮度分量。

在一些示例中，动态成像模块13在获取各个位置的像素序列之前可以对第三干涉图像集中的多张干涉图进行归一化并基于归一化后的干涉图获取各个位置的像素序列。由此，能够对多张干涉图进行规范化。在一些示例中，归一化方法可以包括但不限于是L1归一化或L2归一化。

图6是示出了本公开示例所涉及的带有标定信息D的彩色图像的示意图。在一些示例中，成像系统1可以还包括标定模块(未图示)，标定模块可以用于建立目标频率与颜色的线性对应关系以对目标频率与颜色的对应关系进行标定。在一些示例中，标定模块可以用于获取目标频率的频率范围，并基于频率范围和预设颜色范围之间的线性对应关系对目标频率与颜色的对应关系进行标定。由此，能够直观地获得彩色图像中的颜色对应的频率。在一些示例中，线性对应关系可以表示为：

其中，[H_min，H_max]为频率范围，H_in为输入的频率，H_min为最小的目标频率，H_max为最大的目标频率，C_out为输出的颜色，[C_min，C_max]为预设颜色范围，C_min为预设颜色范围中的最小值，C_max为预设颜色范围中的最大值。在一些示例中，[C_min，C_max]可以为[0，4/3π]。在这种情况下，C_min为0，C_max为4/3π。由此，能够获得更好的视觉效果。但本公开的示例不限于此，在另一些示例中，[C_min，C_max]可以为[0，2π]。作为标定的示例，图6示出了带有标定信息D的鼠视网膜的鲜活样本的彩色图像。如图6所示，标定信息D可以为颜色标尺。在这种情况下，能够直观地获得彩色图像中的颜色对应的频率。

本公开的成像系统1基于奇异值分解获得时间特征向量以通过时间特征向量的波动程度对多张干涉图的伪影进行消除并利用滑动窗对多张干涉图进行信噪比增强处理，最终基于多张干涉图的时域信息和频域信息生成彩色图像。由此，能够降低伪影或者噪声对彩色图像的图像质量的影响，有利于基于彩色图像进行进一步地研究，另外通过本公开生成的彩色图像能够反映生物组织的内部动态信号并且获得该彩色图像的方式是非侵入式。

以下，结合图7详细描述本公开的针对平面干涉仪的信号处理的成像方法。本公开涉及的针对平面干涉仪的信号处理的成像方法有时可以简称为成像方法。成像方法应用于上述的成像系统1中。图7是示出了本公开示例所涉及的针对平面干涉仪的信号处理的成像方法的流程图。

在一些示例中，如图7所示，成像方法可以包括获取多张干涉图并组成第一干涉图像集(步骤S110)，在步骤S110中，可以利用平面干涉仪对样本进行采集以获得多张干涉图。在一些示例中，样本可以是鲜活样本。在一些示例中，多张干涉图可以是通过连续采集获得的。在一些示例中，可以将多张干涉图按照采集时间组成第一干涉图像集。由此，能够方便后续同时对多张干涉图进行处理。具体描述可以参见获取模块10的相关描述。

在一些示例中，如图7所示，成像方法可以包括对第一干涉图像集进行消除伪影处理以获取第二干涉图像集(步骤S120)。在这种情况下，对第一干涉图像集中的伪影进行消除能够提高后续生成的彩色图像的图像质量。在一些示例中，在步骤S120的消除伪影处理中，可以对第一干涉图像集进行奇异值分解以获取时间特征矩阵和奇异值矩阵。在一些示例中，在获取时间特征矩阵和奇异值矩阵后，可以将时间特征矩阵的每列作为时间特征向量。在一些示例中，可以基于时间特征向量的波动程度更新奇异值矩阵中的奇异值。具体地，可以计算时间特征向量的波动程度，基于波动程度筛选符合阈值条件的时间特征向量，然后将该时间特征向量对应的奇异值矩阵中的奇异值更新为预设奇异值。由此，能够基于时间特征向量的波动程度对奇异值矩阵中的奇异值进行更新。在一些示例中，预设奇异值可以为0。具体描述可以参见消除伪影模块11的相关描述。

在一些示例中，在步骤S120中，波动程度可以利用时间特征向量的累积过零率进行表示。也即，可以利用时间特征向量的累积过零率表示时间特征向量的波动程度。如上所述，可以将时间特征矩阵的每列作为时间特征向量。在一些示例中，第i列时间特征向量的累积过零率可以表示为：D_ZRC_i＝|ZRC_i+1-ZRC_i|，其中，i为时间特征矩阵的列索引，ZRC_i为第i列时间特征向量的过零率。由此，能够利用累积过零表示时间特征向量的波动程度。在一些示例中，阈值条件可以为各个时间特征向量的累积过零率大于预设值。在一些示例中，预设值可以为3倍的累积过零率的标准差。由此，能够通过消除累积过零率的离散程度达到预设值程度的信号以消除伪影。具体描述可以参见成像系统1中的消除伪影模块11中关于累积过零率的相关描述。在一些示例中，在步骤S120中，基于更新后的奇异值矩阵可以重建第一干涉图像集并作为第二干涉图像集。具体描述可以参见消除伪影模块11中关于重建第一干涉图像集的相关描述。

在一些示例中，如图7所示，成像方法可以包括对第二干涉图像集进行信噪比增强处理以获取第三干涉图像集(步骤S130)。在步骤S130中，可以利用预设长度的滑动窗对第二干涉图像集中的各个滑动窗中的像素值进行信噪比增强处理以获取第三干涉图像集。在一些示例中，预设长度可以大于1。在一些示例中，预设长度可以为8。在一些示例中，滑动窗的维度与第二干涉图像集的维度可以一致。在一些示例中，滑动窗可以按照预设步长沿着采集时间的方向移动。在一些示例中，预设步长可以为预设长度的预设比例。例如预设比例可以为25％、50％、75％、1/3或2/3等。在一些示例中，在信噪比增强处理中，可以获取各个滑动窗中的像素值的平均值。在一些示例中，可以对各个滑动窗中的像素值与平均值之间的差值求累积和以获取累积值。在一些示例中，可以将累积值的绝对值与预设长度相除以获取平均累积值，将平均累积值作为第三干涉图像集的像素值。由此，能够基于滑动窗和累积和的算法对第二干涉图像集进行信噪比增强处理。具体描述可以参见信噪比增强模块12的相关描述。

在一些示例中，如图7所示，成像方法可以包括基于第三干涉图像集生成彩色图像(步骤S140)。彩色图像可以反映鲜活样本的内部动态信号。在一些示例中，可以基于第三干涉图像集获取三维颜色空间的分量并基于三维颜色空间的分量将第三干涉图像集转为彩色图像。在一些示例中，三维颜色空间的分量可以包括色调分量、饱和度分量和亮度分量(也可以称为明度分量)。在一些示例中，可以获取第三干涉图像集中的多张干涉图的相同位置的多个像素值以形成按照采集时间排序的像素序列，基于各个位置的像素序列的时域信息和频域信息获取各个位置的三维颜色空间的分量，并基于三维颜色空间的分量将第三干涉图像集转为彩色图像。具体描述可以参见动态成像模块13的相关描述。

在一些示例中，在步骤S140中，色调分量可以与目标频率对应。在一些示例中，可以将目标频率作为色调分量。在一些示例中，目标频率可以为平均频率或中值频率。在一些示例中，可以对各个位置的像素序列进行傅立叶变换以获取频率序列，基于频率序列和功率谱密度获取平均频率并作为色调分量。在一些示例中，在步骤S140中，可以基于各个位置的像素序列获取功率谱密度并将功率谱密度的带宽的倒数作为饱和度分量。在一些示例中，可以基于各个位置的像素序列并利用Welch法获取功率谱密度。在一些示例中，在步骤S140中，可以基于各个位置的像素序列获取该像素序列的标准差或方差并将该标准差或该方差作为亮度分量。具体描述可以参见动态成像模块13中关于色调分量、饱和度分量和亮度分量的相关描述。

但本公开的示例不限于此，在另一些示例中，成像方法也可以不包括步骤S120和/或步骤S130。也即，用于生成彩色图像的多张干涉图可以不经过消除伪影处理和/或信噪比增强处理。

在一些示例中，成像方法还包括建立目标频率与颜色的线性对应关系以对目标频率与颜色的对应关系进行标定(未图示)。在一些示例中，可以用于获取目标频率的频率范围，并基于频率范围和预设颜色范围之间的线性对应关系对目标频率与颜色的对应关系进行标定。由此，能够直观地获得彩色图像中的颜色对应的频率。具体描述可以参见标定模块的相关描述。

本公开的成像方法基于奇异值分解获得时间特征向量以通过时间特征向量的波动程度对多张干涉图的伪影进行消除并利用滑动窗对多张干涉图进行信噪比增强处理，最终基于多张干涉图的时域信息和频域信息生成彩色图像。由此，能够降低伪影或者噪声对彩色图像的图像质量的影响，有利于基于彩色图像进行进一步地研究，另外通过本公开生成的彩色图像能够反映生物组织的内部动态信号并且获得该彩色图像的方式是非侵入式。

虽然以上结合附图和实施方式对本发明进行了具体说明，但是可以理解，上述说明不以任何形式限制本发明。本领域技术人员在不偏离本发明的实质精神和范围的情况下可以根据需要对本发明进行变形和变化，这些变形和变化均落入本发明的范围内。

Claims (5)

1.一种干涉图的消除伪影的方法，其特征在于，应用于利用平面干涉仪对鲜活样本的成像中，其中，所述平面干涉仪基于双光束等厚干涉原理进行测量，所述方法包括：获取多张干涉图并按照采集时间组成第一干涉图像集，所述多张干涉图是通过所述平面干涉仪针对所述鲜活样本进行连续采集而获得；对所述第一干涉图像集进行奇异值分解以获取时间特征矩阵和奇异值矩阵，所述时间特征矩阵的每列为时间特征向量；计算所述时间特征向量的波动程度；基于所述时间特征向量的波动程度筛选符合阈值条件的时间特征向量以将该时间特征向量对应的所述奇异值矩阵中的奇异值更新所述奇异值矩阵中的奇异值，并且基于更新后的奇异值矩阵重建所述第一干涉图像集并作为第二干涉图像集，利用所述时间特征向量的累积过零率表示所述时间特征向量的波动程度，第i列时间特征向量的累积过零率表示为：D_ZRC_i＝|ZRC_i+1-ZRC_i|，其中，i为所述时间特征矩阵的列索引，ZRC_i为第i列时间特征向量的过零率；利用预设长度的滑动窗对所述第二干涉图像集中的各个滑动窗中的像素值进行信噪比增强处理以获取第三干涉图像集，所述滑动窗按照预设步长沿着所述采集时间的方向移动，所述滑动窗的维度与所述第二干涉图像集的维度一致；基于所述第三干涉图像集生成彩色图像，所述彩色图像反映所述鲜活样本的内部动态信号并且获得所述彩色图像的方式是非侵入式，在基于所述第三干涉图像集生成所述彩色图像中，基于所述第三干涉图像集获取三维颜色空间的分量并基于所述三维颜色空间的分量将所述第三干涉图像集转为所述彩色图像，在基于所述第三干涉图像集获取所述三维颜色空间的分量中，获取所述第三干涉图像集中的多张干涉图的相同位置的多个像素值以形成按照所述采集时间排序的像素序列，基于各个位置的像素序列的时域信息和频域信息获取各个位置的所述三维颜色空间的分量，所述三维颜色空间的分量包括色调分量、饱和度分量和亮度分量，基于所述频域信息获取目标频率，所述色调分量与所述目标频率对应。

2.根据权利要求1所述的方法，其特征在于：

所述阈值条件为各个时间特征向量的累积过零率大于预设值，其中，所述预设值为3倍的所述累积过零率的标准差。

3.根据权利要求1所述的方法，其特征在于：

绘制各个时间特征向量的曲线并基于所述曲线求该时间特征向量的过零率。

4.根据权利要求1所述的方法，其特征在于：

对转成二维数据的所述第一干涉图像集进行所述奇异值分解以获取所述时间特征矩阵和所述奇异值矩阵，所述二维数据的大小为m×n，m为一张干涉图的尺寸，n为所述多张干涉图的数量。

5.根据权利要求1所述的方法，其特征在于：

在所述信噪比增强处理中，获取各个滑动窗中的像素值的平均值，对各个滑动窗中的像素值与所述平均值之间的差值求累积和以获取累积值，将所述累积值的绝对值与所述预设长度相除以获取平均累积值，将所述平均累积值作为所述第三干涉图像集的像素值。