CN105245761A - 单次快照多频解调方法 - Google Patents

Abstract

本发明涉及一种单次快照多频解调的方法，在时间域或空间域对由一个或多个原始成分经不同频率调制相加得到的调制图像、尤其是包含有多个不同频率成分的调制图像，首先依次提取出每个像素点的各个频率交流(AC)和直流(DC)分量值，进而得到相应各个频率交流(AC)和直流(DC)的原始成分图像。该方法可以用于时域或空域中，利用一次测量就能分解出多个频率成分图像，具有速度快、解调精度更高以及去噪效果很好的优点，满足了一次性获取多个频率信息量的需求，克服了多次测量中不能避免的误差问题。同时，也可以利用该解调方法一次性传输多幅图像信息，实现在通讯领域多信息并行实时传输。

Description

单次快照多频解调方法

技术领域

本发明涉及空间频域成像、实时多成分成像、时空域实时信号调制解调技术及多组分图像信息传输，更具体地，涉及一种用于非接触式成像技术中涉及的单次快照多频解调方法。

背景技术

在生物医学成像领域中，新兴的空间频域成像(SpatialFrequencyDomainImaging，SFDI)作为一种新颖的非接触式成像技术，具有独一无二的空间上同时解析光学吸收和散射参数的能力，允许宽视场量化组织光学参数分布。通过入射不同空间频率的空间调制图案到样品区域并用CCD相机捕捉反射图像来获取样品的调制传递函数(MTF)。MTF包含了重要的光学特性信息——吸收系数(μ_a)和衰减散射系数(μ_s')。基于蒙特卡罗或各种散射模型，通过非线性最小二乘法拟合或查表方法，由MTF数据可以反演计算出生物组织的吸收系数与衰减散射系数的二维分布图。最后，通过光学参数的变化可以反推出组织结构和组织成分的改变，进而做出相应疾病的诊断。

根据EssexT.J.H.，ByrneRO.AlaserDopplerscannerforimagingbloodflowinskin[J].Medicalengineeringandphysics，1991，13(3)：189-194的记载，假设入射到样品的结构光强用函数表示为：

$S=\frac{S_0}{2}\[1+M_{0}cos(2{\mathrm{\πf}}_{x}x+\mathrm{\α})\]---\left(1\right)$

这里S₀代表光源强度，M₀是入射调制深度，f_x是空间频率，α是空间相位，x是空间坐标。

从样品反射出来被CCD照相机捕获的光可以分解为直流DC部分和交流AC部分：

I＝I_AC+I_DC(2)

样品反射光的AC部分可以用函数表示为：

I_AC＝M_AC(x,f_x)×cos(2πx+α)(3)

这里M_AC表征散射光子密度波的调制，此因子在混沌介质中依赖于组织的光学特性，目前的主流作法是基于扩散理论或蒙特卡罗的光传输方法进行建模。为了得出M_AC，必须对信号进行解调，传统的标准作法是三相移法(在文献NeilMAA,JuskaitisR,WilsonT.Methodofobtainingopticalsectioningbyusingstructuredlightinaconventionalmicroscope.Opt.Lett1997；22(24):1905–1907.[PubMed:18188403]中提到)。即若样品是在一个特定频率正弦波的三个相位差α＝0,2π/3,4π/3处进行光照，测得三幅光强图像I₁、I₂、I₃，则M_AC因子可以用式(4)解调方程来计算。

$M_{AC}(x,f_x)=\frac{\sqrt{2}}{3}{\[{(I_1-I_2)}^2+{(I_2-I_3)}^2+{(I_3-I_1)}^2\]}^{\frac{1}{2}}---\left(2\right)$

为了进一步得到组织的光学参数，M_AC需要在不同空间频率使用三个相位投影来测量。首先，多种频率的光投射到样品，这种方式是以多个相位投影到样品并用方程(4)进行解调。然后，漫反射在每一个空间频率都使用硅校准模型的已知光学参数来进行校准，从而矫正MTF值。最后在图像上的每一个像素点使用一个逆模型得出各独立波长的光学参数。

总的来说，传统的SFDI的成像步骤及光学参数的获得步骤如下：

a)将包含多种频率f_x的调制光投影到样品上，并通过CCD相机采集从样品反射回的光强；

b)每一个光照频率在三个相位点成像，然后用解调公式(4)进行解调，并由公式(5)得出每一个像素点的反射率R，其中MTF_system由已知光学参数校准模型在同一条件下测量获得：

M_AC(x_i)＝I₀MTF_system(x_i)×R(x_i)(5)

c)每一个像素点的R值利用光传输模型的蒙特卡罗方法或查表获得吸收系数μ_a和衰减散射系数μ_s′的二维映射分布。

由上可知，三相移标准法，即给定三个不同初始相位(0°,120°,240°)，通过公式求解出交流分量和直流分量，该方法被公认为解调交流/直流分量的“金标准”，但该方法在实际成像中需要至少三次成像才能够解调出交流分量，限制了成像时间和成像次数。除此之外，根据Nadeau，K.P.，Durkin，A.J.，Tromberg，B.J..Advanceddemodulationtechniquefortheextractionoftissueopticalpropertiesandstructuralorientationcontrastinthespatialfrequencydomain[J].JournalofBiomedicalOptics，2014，19(5)：056013.的记载，还可以利用Hilbert变换法在单个相位的条件下解调出交流分量，这能够大大提高光学参数测量效率，但其仅能实现单个相位的交流分量解调，并且噪声抑制效果不佳。

发明内容

本发明的目的在于克服现有技术存在的上述缺陷，提供了一种单次快照多频解调法，该方法可以快速解调出多个不同频率的交流成分幅值和直流分量。

为达到上述目的，本发明采取的解决方案为：

一种单次多频快照解调法，其包括：在时间域或空间域中，对于包含有一个或多个频率的调制图像，首先依次提取出每个像素点的各个频率交流(AC)和直流(DC)分量值，进而得到相应于各个频率AC和DC的原始成分图像。

进一步地，对于具有一种频率或两种以上不同频率的交流分量的调制图像，首先，取一T₁×T₂的矩阵作为基核(Kernel)，其中T₁、T₂分别为各交流分量周期的横向和纵向分量最小公倍数；然后在图像中，按照预定顺序依次以每个像素点为矩阵左上角选取基核大小的基核单元，并分别与相对应的基核图案作积分求和处理。

进一步地，其中积分求和处理进一步包括：图像中的每个基核单元分别与同频率、同方向的余弦及正弦基核图案作乘法处理后，再进行积分求和，从而得到与原始图像尺寸大小相同的余弦调制和正弦调制图案。

进一步地，其中对于具有不同频率f₁、f₂、…的交流分量的调制图像，各分量关系如公式(6)所示：

$\sqrt{m_i^2+n_i^2}=1---\left(7\right)$

其中，k≥1为交流调制分量个数，A_i、f_i和分别为各交流成分幅值、频率和空间初始相位，m_i∈[-1，1]、n_i∈[-1,1]共同决定了交流调制条纹的方向，且需满足公式(7)，B为直流分量，x、y为空间坐标。

其中，进一步利用公式(8)解调出调制图案的交流成分幅值，利用公式(9)解调出直流成分：

$A_i=\frac{\sqrt{{\[\underset{\mathrm{\σ}}{\∫\∫}f(x,y)\·\cos (2{\mathrm{\πf}}_i\·(m_i\·x+n_i\·y\left)\right)dxdy\]}^2+{\[\underset{\mathrm{\σ}}{\∫\∫}f(x,y)\·\sin (2{\mathrm{\πf}}_i\·(m_i\·x+n_i\·y\left)\right)dxdy\]}^2}}{{\[\underset{\mathrm{\σ}}{\∫\∫}2{\mathrm{cos\πf}}_i\·(m_i\·x+n_i\·y)dxdy\]}^2}---\left(8\right)$

$B=\frac{1}{T_1\×T_2}\underset{\mathrm{\σ}}{\∫\∫}f(x,y)dxdy---\left(9\right)$

其中，σ为图像中每个像素点对应的基核单元。T₁、T₂分别为各交流分量周期的横向和纵向分量最小公倍数；

其中，进一步地，对于基核图案，提取图像的不同频率或方向的交流成分时，基核图案也不相同，并且基核图案与提取交流调制图案具有相同的频率和方向。

其中，进一步地，所述预定顺序为以图像左上角为起点，依次从左到右、从上到下的顺序。

此外，还提供了一种空间频域成像方法，其包括将包含一种或多种频率的调制光投影到样品上，并通过CCD相机采集从样品反射回的光强、透射出光强或激发出的荧光光强；其中所述CCD相机采集的光分解为直流(DC)部分和交流(AC)部分，其特征在于所述CCD相机采集的光采用前面所述的单次多频快照解调法进行解调。

此外还提供了一种信号传输方法，其中包括前述的空间频域信息的调制及解调方法。

由此可见，本发明提出一种单次快照多频解调法，在空间域由一个或多个原始成分经不同频率调制相加得到的调制图像、尤其是包含有多个不同频率成分的调制图像，通过该方法依次提取出每个像素点的各个频率交流(AC)和直流(DC)分量值，进而得到相应于各个频率交流(AC)和直流(DC)的原始成分图像。该方法具有速度快、解调精度更高以及去噪效果很好的优点，满足了一次性获取多个频率信息量的需求，克服了多次测量中不能避免的误差问题。同时，也可以利用该解调方法一次性传输多幅图像信息，实现在通讯领域多信息并行实时传输。

附图说明

为了更清楚地说明本发明实施例的技术方案，下面将对实施例或现有技术描述中所需要使用的附图作简单地介绍，显而易见地，下面描述中的附图仅仅是本发明的一些示例性实施例。

图1为利用单次快照多频解调法对包括两个空间频率的图像进行解调处理的图示。

图2为分别利用三相移标准法、Hilbert变换法和单次快照多频解调法对含有一个频率的图案进行解调处理的对比图示。

图3为分别利用Hilbert变化法和单次快照分解法对含有两个频率的图案进行解调处理的对比图示。

图4为利用单次快照分解法对含有三个频率的图案进行解调处理的对比图示。

具体实施方式

下面将结合本发明的附图，对本发明实施方式中的技术方案进行清楚、完整地描述，显然，所描述的实施方式仅仅是本发明的一部分实施方式，而不是全部的实施方式。基于本发明中的实施方式，本领域普通技术人员在没有作出创造性劳动前提下所获得的所有其他实施方式，都属于本发明保护的范围。

下面就以空间频域成像为例介绍了单次快照多频解调法，相较于标准三相移法，速度快、解调精度更高并具有很好的去噪效果，满足了一次性获取多个频率信息量的需求，克服了多次测量中不能避免的误差问题。

单次多频快照解调法是在空间域对由一个或多个原始成分经不同频率调制相加得到的调制图像通过该方法依次提取出每个像素点的各个频率AC和DC分量值，进而得到相应于各个频率AC和DC的原始成分图像。如图1所示，利用单次快照多频解调法对包括两个空间频率的图像进行解调处理，单次测得的图像中包括两个不同空间频率(f₁、f₂)的成分。单次快照多频解调法可以精准地解调出不同频率的空间调制图像。

具体步骤如下：

对于具有两种不同频率交流分量的图像，若频率分别为f₁、f₂，各分量关系如公式(10)所示。为了能够提取出图像中不同频率交流分量大小，首先，取一T₁×T₂的矩阵作为基核(Kernel)，其中T₁、T₂分别为各交流分量周期的横向和纵向分量最小公倍数；然后在图像中，依次从左到右、从上到下的顺序以每个像素点为矩阵左上角选取基核大小的基核单元，并分别与相对应的基核图案作积分求和处理。对于基核图案，提取图像的不同频率、方向交流分量时，基核图案也不相同，但应保证基核图案与提取分量具有相同的频率和方向。图像中的每个基核单元分别与同频率、同方向的余弦及正弦图案作乘法处理后，进行积分求和。得到与原始图案尺寸大小相同的余弦解调和正弦解调图案；

然后，利用公式(11)、(12)解调出调制图案的交流成分幅值，利用公式(13)解调出直流成分。

其中在公式(10)中，a、b分别为各交流成分幅值，f₁和f₂为各交流调制图案频率，和为空间初始相位，c为直流分量，x、y为空间坐标。

$a=\frac{\sqrt{{\[\underset{\mathrm{\σ}}{\∫\∫}f(x,y)\cos (2{\mathrm{\πf}}_1\·x)dxdy\]}^2+{\[\underset{\mathrm{\σ}}{\∫\∫}f(x,y)\sin (2{\mathrm{\πf}}_1\·x)dxdy\]}^2}}{\[\underset{\mathrm{\σ}}{\∫\∫}\cos (2{\mathrm{\πf}}_1\·x)dxxy\]}---\left(11\right)$

$b=\frac{\sqrt{{\[\underset{\mathrm{\σ}}{\∫\∫}f(x,y)cos(2{\mathrm{\πf}}_2\·(x+y\left)\right)dxdy\]}^2+{\[\underset{\mathrm{\σ}}{\∫\∫}f(x,y)sin(2{\mathrm{\πf}}_2\·(x+y\left)\right)dxdy\]}^2}}{{\[\underset{\mathrm{\σ}}{\∫\∫}\cos (2{\mathrm{\πf}}_2\·(x+y\left)\right)dxdy\]}^2}$

$\left(12\right)$

$c=\frac{1}{T_1\×T_2}\underset{\mathrm{\σ}}{\∫\∫}f(x,y)dxdy---\left(13\right)$

其中，σ为图像中每个像素点对应的基核单元。T₁、T₂分别为各交流分量周期的横向和纵向分量最小公倍数。

进一步地，对于具有多个不同频率f₁、f₂、…的交流分量的调制图像，各分量关系则如公式(14)所示：

其中，k≥1为交流调制分量个数，A_i、f_i和分别为各交流成分幅值、频率和空间初始相位，m_i∈[-1，1]、n_i∈[-1,1]共同决定了交流调制条纹的方向，B为直流分量，x、y为空间坐标；并且其中调制图案的交流成分幅值利用公式(15)解调出，直流成分利用公式(16)解调出：

$A_i=\frac{\sqrt{{\[\underset{\mathrm{\σ}}{\∫\∫}f(x,y)\·\cos (2{\mathrm{\πf}}_i\·(m_i\·x+n_i\·y\left)\right)dxdy\]}^2+{\[\underset{\mathrm{\σ}}{\∫\∫}f(x,y)\·\sin (2{\mathrm{\πf}}_i\·(m_i\·x+n_i\·y\left)\right)dxdy\]}^2}}{{\[\underset{\mathrm{\σ}}{\∫\∫}\cos (2{\mathrm{\πf}}_i\·(m_i\·x+n_i\·y\left)\right)dxdy\]}^2}---\left(15\right)$

$B=\frac{1}{T_1\×T_2}\underset{\mathrm{\σ}}{\∫\∫}f(x,y)dxdy---\left(16\right)$

其中，σ为图像中每个像素点对应的基核单元。T₁、T₂分别为各交流分量周期的横向和纵向分量最小公倍数。

对比实验：

分别利用三相移标准法、Hilbert变换法和单次快照多频解调法对含有一个频率的图案进行解调，得到的结果如图2所示。其中原始图像包含三部分分量，分别为直流分量、交流分量及噪声信号，我们通过几种不同的解调方法来比较各方法的优缺点。根据所得结果可知，三种方法都可实现从噪声及直流信号中解调出交流信号的目的，但从对噪声抑制效果来说，单次快照多频解调法优于三相移标准法和Hilbert变换法，其能够很好地抑制噪声干扰。

如图3所示，利用Hilbert变化法和单次快照分解法分别对含有两个频率的图案进行解调：其中图3(a).原始图像是通过对“中”字作横向调制，f_y＝0.2，对“华”字作纵向调制，f_x＝0.1；然后对“中”和“华”字进行叠加处理，并加以直流分量和噪声信号。已知三相移标准法只能解调出单频交流分量，因此在此只对Hilbert变换法和单次快照多频解调法作比较。图3(b)为对图3(a)作Hilbert变换法解调后的结果，从结果可以看出，在有两个或多个频率交流分量情况下，解调后的交流分量并不能分离出不同频率的信号，而从图3(c)与3(d)可以看出，单次快照多频解调法能够很好地解调出不同频率的交流调制信号。

如图4所示，利用单次快照分解法分别对含有三个频率的图案进行解调：其中图4(a)中，原始图像是通过对“大”字作横向调制，f_y＝0.1，对“中”字作纵向调制，f_x＝0.5，对“华”作倾斜45°调制，f_x＝0.25，f_y＝0.25；然后对“大”、“中”和“华”字进行叠加处理，并加以直流成分和噪声信号。已知三相移标准法只能解调出单频交流分量，Hilbert变换法区分不出多个频率调制信号，下图仅为单次快照多频解调法得到的结果。从图4(b)、图4(c)与图4(d)可以看出，单次快照多频解调法能够很好地解调出不同频率的交流调制信号。

从实验结果可以清楚得到，单次快照多频解调法不仅可以很好地从噪声中解调出单频交流分量及抑制噪声的干扰，而且可以解调出两个或多个频率的交流分量。这样我们能够实现不同空间频率MTF函数的实时成像，快速解析出组织的光学参数，从而克服多次成像引起噪声干扰的难题。

如上所述的各实施方式，通过利用上述单次快照多频解调法，在空间域对由一个或多个原始成分经不同频率调制相加得到的调制图像、尤其是包含有多个不同频率成分的调制图像，通过该方法依次提取出每个像素点的各个频率交流(AC)和直流(DC)分量值，进而得到相应于各个频率AC和DC的原始成分图像。该方法具有速度快、解调精度更高以及去噪效果很好的优点，满足了一次性获取多个频率信息量的需求，克服了多次测量中不能避免的误差问题。同时，也可以利用该解调方法一次性传输多幅图像信息，实现在通讯领域多信息并行实时传输。

本技术领域技术人员可以理解，本发明中已经讨论过的各种操作、方法、流程中的步骤、措施、方案可以被交替、更改、组合或省略。进一步地，具有本发明中已经讨论过的各种操作、方法、流程中的其他步骤、措施、方案也可以被交替、更改、重排、分解、组合或省略。进一步地，现有技术中的具有与本发明中公开的各种操作、方法、流程中的步骤、措施、方案也可以被交替、更改、重排、分解、组合或省略。

以上所述仅是本发明的部分实施方式，应当指出，对于本技术领域的普通技术人员来说，在不脱离本发明原理的前提下，还可以做出若干改进和调整，这些改进和调整也落入本发明的保护范围内。

Claims (9)

1.一种单次多频快照解调方法，其包括：

在时间域或空间域中，对于包含有一个或多个频率的调制图像，首先依次提取出每个像素点的各个频率交流(AC)和直流(DC)分量值，进而得到相应于各个频率AC和DC的原始成分图像。

2.根据权利要求1所述的方法，具体包括：对于所述调制图像，首先，取一T₁×T₂的矩阵作为基核(Kernel)，其中T₁、T₂分别为各交流分量周期的横向和纵向分量最小公倍数；然后在图像中，按照预定顺序依次以每个像素点为矩阵左上角选取基核大小的基核单元，并分别与相对应的基核图案作积分求和处理。

3.根据权利要求2所述的方法，其中积分求和处理进一步包括：图像中的每个基核单元分别与同频率、同方向的余弦及正弦基核图案作乘法处理后，再进行积分求和，从而得到与原始图像尺寸大小相同的余弦解调和正弦解调图案。

4.根据权利要求1至3任一项所述的方法，其中对于具有各不同频率f₁、f₂、…的交流分量的调制图像，各分量关系如公式(1)所示：

$\sqrt{m_i^2+n_i^2}=1---\left(2\right)$

其中，k≥1为交流调制分量个数，A_i、f_i和分别为各交流成分幅值、频率和空间初始相位，m_i∈[-1，1]、n_i∈[-1,1]并共同决定了交流调制条纹的方向，且需满足公式(2)，B为直流分量，x、y为空间坐标。

5.根据权利要求4所述的方法，其中进一步利用公式(3)解调出调制图案的交流成分幅值，利用公式(4)解调出直流成分：

$A_i=\frac{\sqrt{{\[\underset{\mathrm{\σ}}{\∫\∫}f(x,y)\·\cos (2{\mathrm{\πf}}_i\·(m_i\·x+n_i\·y\left)\right)dxdy\]}^2+{\[\underset{\mathrm{\σ}}{\∫\∫}f(x,y)\·\sin (2{\mathrm{\πf}}_i\·(m_i\·x+n_i\·y\left)\right)dxdy\]}^2}}{{\[\underset{\mathrm{\σ}}{\∫\∫}\cos (2{\mathrm{\πf}}_i\·(m_i\·x+n_i\·y\left)\right)dxdy\]}^2}---\left(3\right)$

$B=\frac{1}{T_1\×T_2}\underset{\mathrm{\σ}}{\∫\∫}f(x,y)dxdy---\left(4\right)$

其中，σ为图像中每个像素点对应的基核单元，T₁、T₂分别为各交流分量周期的横向和纵向分量最小公倍数。

6.根据权利要求2至5任一所述的方法，其中对于基核图案，提取图像的不同频率或方向的交流成分时，基核图案也不相同，并且基核图案与提取交流调制图案具有相同的频率和方向。

7.根据权利要求2至6任一所述的方法，其中所述预定顺序为以图像左上角为起点，依次从左到右、从上到下的顺序。

8.一种空间频域成像方法，其包括将包含一种或多种频率的调制光投影到样品上，并通过CCD相机采集从样品反射回的光强、透射出光强或激发出的荧光光强；其中所述CCD相机采集的光分解为直流(DC)部分和交流(AC)部分，其特征在于所述CCD相机采集的光采用权利要求1至7任一项所述的单次多频快照解调法进行解调。

9.一种信号传输方法，其中包括权利要求1-8所述的空间频域信息的调制及解调方法。