CN101966078B - 一种近红外漫射光频域信息获取方法 - Google Patents

Abstract

本发明属于组织光学研究中的光学参数测量领域，涉及一种近红外漫射光频域信息获取方法，包括下列步骤：首先建立一个双距离探测频域系统；利用频域蒙特卡洛模拟方法，计算系统频域相位延迟参数Φ_cal：令受检体为标准体，改变探测条件，得到四个标准体的输出相位；令受检体为未知组织体，与步骤(3)的探测过程相同，得到四对未知组织体的输出幅度和输出相位；计算未知组织体在频域中的幅值信息；计算未知组织体在频域中的相位信息。本发明的获取近红外漫射光频域信息方法，能够消除系统固有幅度衰减、固有相位延迟和测试操作过程方法等因素影响。

Description

一种近红外漫射光频域信息获取方法

技术领域

本发明属于组织光学中的光学参数测量领域，具体涉及一种用于组织光学研究的近红外漫射光频域信息获取方法。

背景技术

频域技术⁽¹⁾是近红外漫射光检测技术之一，因其具有可测信息丰度大、测量精度较高、设备价格较低的优点而获得研究者的广泛关注。所谓近红外漫射光检测技术，即利用特定波长(600-900nm)的光照射组织体，检测到达组织体表以下数厘米后出射光的分布情况，并根据检测结果重构出被检组织的光学特性参数(吸收系数μ_a，约化散射系数μ′_s)，从而达到组织病变的分级以及监测癌变组织对治疗的反应的目的；一般而言，癌变组织的吸收系数μ_a比正常组织高，约化散射系数μ′_s比正常组织低。组织体吸收系数和约化散射系数的变化表现为，同样的入射光源情况下出射光的状态不同(频域技术中表现为测得的幅度衰减和相位延迟不同)。在频域技术中，先利用频域检测系统测得光强被高频信号调制的入射光源经过待测组织体后的交流幅度衰减AC和相位延迟Φ，再利用AC和Φ重构出组织体的光学参数即吸收参数μ_a及漫反射参数μ′_s，最后根据μ_a和μ′_s达到确诊所测组织体的生理病理状态的目的。由此可见，对组织性质的如实诊断和分辨最终归结到幅度AC和相位Φ的精确测量，其过程需要有效的系统校正方法的保证。

近红外漫射光频域检测系统校正的目的是消除由测量系统中固有的光路元器件、电路元器件、测量操作过程等因素产生的幅度衰减和相位延迟，得到仅由被检组织体引起的幅度和相位变化。通常采用直接绝对校正法⁽²⁾、空探测参考路校正法^(3，4)以及光学参数已知的标准体校正法⁽⁵⁾，其中空探测参考路校正法包含直接绝对校正法的思想。直接绝对校正法系统配置简单；空探测参考路校正法可以有效降低由光源引起的系统漂移和高频电路噪声，但无法消除参考通道和测量通道之间相位和幅度的串扰，并且对参考路和测试路的测量对称性要求较高，理论上要求参考路探测器和测量路探测器性能完全一致；光学参数已知的标准体校正法要求在进行标准体和未知组织体测量时入射光源功率相同，因此实际选择标准体时需使标准体的光学参数与被测组织相近，以保证在相同的输入光功率下标准体测量和未知组织体测量的数据结果均能获得足够的信噪比，还要求标准体的几何结构与待测组织体一致以降低测试操作过程对数据结果的影响。对于组织体漫反射光的测量，上述三种方法都无法有效克服测量过程中光纤与被测组织之间耦合因素的影响。

参考文献：

(1)Sergio Fantini，Beniamino Barbieri，Maria Angela Franceschini，Enrico Gratton，Frequency-domainspectroscopy，Applications of Optical Engineering to the Study of Cellular Pathology，Vol.1，p.57-66，1997.

(2)Elizabeth M.C.Hillman，Jeremy C.Hebden，Florian E.W.Schmidt，Simon R.Arridge，MartinSchweiger，Hamid Dehghani，David T.Delpy，Calibration techniques and datatype extraction for time-resolvedoptical tomography，Rev.Sci.Instrum.Vol.71(9)，p.3415-3427，2000.

(3)Ilkka Nissila¨，I Kalle Kotilahti，Kim Fallstro¨m，Toivo Katila，Instrumentation for the accuratemeasurement of phase and amplitude in optical tomography，Rev.Sci.Instrum.Vol.73(9)，p.3306-3312，2002.

(4)Joshua B.Fishkin，Peter T.C.So，Albert E.Cerussi，Sergio Fantini，Maria Angela Franceschini，EnricoGratton，Frequency-domain method for measuring spectral properties in multiple-scattering media：methemoglobin absorption spectrum in a tissuelike phantom，APPLIED OPTICS，Vol.34(7)，p.1143-1155，1995.

(5)Sonja Spichtig，RenéHornung，Derek W.Brown，Daniel Haensse，Martin Wolf，Multifrequencyfrequency-domain spectrometer for tissue analysis，Rev.Sci.Instrum.Vol.80(024301)，2009.

发明内容

本发明的目的是，克服现有技术的上述不足，提出一种近红外漫射光频域信息获取方法，该方法采用一种标准体-双频率-双距离校正方法，能够消除系统固有幅度衰减、固有相位延迟和测试操作过程方法等因素影响，这里的“双频率”中频率表示光源调制频率，“双距离”在距离表示入射光点与探测点之间的距离。本发明采用如下的技术方案：

一种近红外漫射光频域信息获取方法，包括下列步骤：

(1)首先建立一个双距离探测频域系统，在该系统中，光源为在频率可调的近红外光源，用于为系统提供入射到被检组织的光源，光源发出的光分成两束，分别被导入到受检体的两个不同位置，被受检体在探测点处的漫反射光由探测器检测出，设两束光入射到受检体的两个不同位置距离探测点的距离分别为探测距离d₁和探测距离d₂，设系统的输入幅度为AC_in，初始相位为Φ_in，输出幅度为AC_out，输出相位为Φ_out；

(2)利用频域蒙特卡洛模拟得到光学参数已知的标准体在光源调制频率为f₁、探测距离为d₁时的理论相位以及探测距离为d₂时的理论相位

利用频域蒙特卡洛模拟得到光学参数已知的标准体在光源调制频率为f₂、探测距离为d₁时的理论相位

以及探测距离为d₂时的理论相位

并根据下列公式计算系统频域相位延迟参数Φ_cal：

(3)令受检体为标准体，分别探测光源调制频率为f₁、探测距离为d₁时标准体的输出相位

光源调制频率为f₁、探测距离为d₂时标准体的输出相位

探测光源调制频率为f₂、探测距离为d₁时标准体的输出相位

和探测光源调制频率为f₂、探测距离为d₂时的标准体输出相位

(4)令受检体为未知组织体，分别探测光源调制频率为f₁、探测距离为d₁时的输出幅度AC_out(f₁，d₁)，光源调制频率为f₁、探测距离为d₂时的输出幅度

和输出相位

探测光源调制频率为f₂、探测距离为d₁时的输出幅度

和输出相位

以及探测光源调制频率为f₂、探测距离为d₂时的输出幅度

(f₂，d₂)和输出相位

(5)根据下列公式计算未知组织体在频域中的幅值信息A^组织体：

(6)根据下列公式计算未知组织体在频域中的相位信息Θ^组织体：

本发明利用光学参数已知的标准体、两个光源强度调制频率、两个探测点进行频域系统幅度和相位的测试；设计合理的测试流程及测试数据处理方法，消除系统固有因素及测量方法过程(如光纤与被测体间耦合因素)对数据结果的影响，得到仅由受检体衰减后的幅度和仅由受检体延迟后的相位。达到有效校正系统提高频域系统参数测量准确性、最终提高被检组织光学参数重构准确率的目的。具体而言，本发明的频域检测方法的主要特点体现在：

1、本发明仅相位校正时需参考标准体的测量数据，幅度校准不需要标准体的测量数据；

2、本发明校正过程中不需要测量标准体时的入射光功率和测量待检组织体时的入射光功率相同，不需要标准体和待检组织体的光学参数相近；

3、本发明可以消除光纤与光开关、组织体等之间接触耦合因素的影响，不要求标准体和待检组织体的几何结构相似。

附图说明

图1：本发明适用的频域检测系统的测量原理简图。

图2：系统幅度衰减过程示意图。

图3：系统相位延迟过程示意图。

具体实施方式

本发明通过对标准体进行双频率双距离测试，能够得到系统固有的相位延迟和幅度衰减；在此基础上测量待检组织体能够得到仅由组织体引起的幅度和相位变化，提高频域系统参数测量的准确性，从而能够提高组织体光学参数重构的精确性。下面根据各个附图和实施例对本发明做进一步描述。

根据频域系统原理，本发明适用的频域检测系统的测量原理简图如图1所示：

(1)近红外波段范围内的光源1，用于为系统提供入射到被检组织的光源；

(2)双通道光开关2，用于将由光纤[1]传导的入射光源通过光纤[2-1]和光纤[2-2]分时导入到受检体(标准体或组织体)的两个不同位置，使双距离(d₁，d₂)探测得以实现；

(3)探测器3，用于检测探测点处由光纤[3]传导的反映受检体吸收、散射性质的漫反射光。

图1中引起交流幅度衰减的过程如图2，包括：

(1)η₁表示光纤[1]产生的幅度衰减；

(2)η₂表示光纤[2-1]与光开关接触耦合时产生的幅度衰减，η₃表示光纤[2-1]自身产生的幅度衰减，η₄表示光纤[2-1]与受检体耦合接触时产生的幅度衰减；

(3)η′₂表示光纤[2-2]与光开关接触耦合时产生的幅度衰减，η′₃表示光纤[2-2]自身产生的幅度衰减，η′₄表示光纤[2-2]与受检体耦合接触时产生的幅度衰减；

(4)X(d₁)表示探测距离为d₁时受检体产生的幅度衰减，X(d₂)表示探测距离为d₂时受检体产生的幅度衰减；

(5)η₅表示光纤[3]与受检体接触耦合时产生的幅度衰减，η₆表示光纤[3]自身产生的幅度衰减。

图1中引起相位延迟的过程如图3，包括：

(1)θ₁表示光纤[1]产生的相位延迟；

(2)θ₂表示光纤[2-1]与光开关接触耦合所产生的相位延迟，θ₃表示光纤[2-1]自身产生的相位延迟，θ₄表示光纤[2-1]与受检体接触耦合所产生的相位延迟；

(3)θ′₂表示光纤[2-2]与光开关接触耦合所产生的相位延迟，θ′₃表示光纤[2-2]自身产生的相位延迟，θ′₄表示光纤[2-2]与受检体接触耦合所产生的相位延迟；

(4)Φ(d₁)表示探测距离为d₁时受检体产生的相位延迟，Φ(d₂)表示探测距离为d₂时受检体产生的相位延迟；

(5)θ₅表示光纤[3]与受检体接触耦合所产生的相位延迟，θ₆表示光纤[3]自身产生的相位延迟。

设系统的输入幅度为AC_in，初始相位为Φ_in，输出幅度为AC_out，输出相位为Φ_out；以下推导过程假设实际中图2中的幅度衰减和图3中的由光纤耦合因素引起的相位延迟均与光源调制频率无关。

本发明的频域参数获取的理论依据如下：

1、设光源调制频率为f₁，探测距离为d₁，探测得到，

AC_out(f₁，d₁)＝AC_in(f₁)·η₁·η₂·η₃·η₄·X(f₁，d₁)·η₅·η₆             (1)

Φ_out(f₁，d₁)＝Φ_in(f₁)-θ₁(f₁)-θ₂-θ₃(f₁)-θ₄-Φ(f₁，d₁)-θ₅-θ₆(f₁)       (2)

2、保持光源调制频率不变，改变探测距离为d₂，探测得到，

AC_out(f₁，d₂)＝AC_in(f₁)·η₁·η′₂·η′₃·η′₄·X(f₁，d₂)·η₅·η₆       (3)

Φ_out(f₁，d₂)＝Φ_in(f₁)-θ₁(f₁)-θ′₂-θ′₃(f₁)-θ′₄-Φ(f₁，d₂)-θ₅-θ₆(f₁) (4)

由(1)、(3)得，

$\frac{{\mathrm{AC}}_{\mathrm{out}}(f_1,d_2)}{{\mathrm{AC}}_{\mathrm{out}}(f_1,d_1)}=\frac{{\mathrm{\η}}_2^{\′}\·{\mathrm{\η}}_3^{\′}\·{\mathrm{\η}}_4^{\′}\·X(f_1,d_2)}{{\mathrm{\η}}_2\·{\mathrm{\η}}_3\·{\mathrm{\η}}_4\·X(f_1,d_1)}---\left(5\right)$

由(2)、(4)得，

Φ_out(f₁，d₂)-Φ_out(f₁，d₁)＝θ₂+θ₃(f₁)+θ₄+Φ(f₁，d₁)-θ′₂-θ′₃(f₁)-θ′₄-Φ(f₁，d₂)  (6)

3、将光源调制频率改为f₂，调整探测距离为d₁，探测得到，

AC_out(f₂，d₁)＝AC_in(f₂)·η₁·η₂·η₃·η₄·X(f₂，d₁)·η₅·η₆              (7)

Φ_out(f₂，d₁)＝Φ_in(f₂)-θ₁(f₂)-θ₂-θ₃(f₂)-θ₄-Φ(f₂，d₁)-θ₅-θ₆(f₂)        (8)

4、保持光源状态f₂不变，改变探测距离为d₂，探测得到，

AC_out(f₂，d₂)＝AC_in(f₂)·η₁·η′₂·η′₃·η′₄·X(f₂，d₂)·η₅·η₆        (9)

Φ_out(f₂，d₂)＝Φ_in(f₂)-θ₁(f₂)-θ′₂-θ′₃(f₂)-θ′₄-Φ(f₂，d₂)-θ₅-θ₆(f₂)  (10)

由(7)，(9)得

$\frac{{\mathrm{AC}}_{\mathrm{out}}(f_2,d_2)}{{\mathrm{AC}}_{\mathrm{out}}(f_2,d_1)}=\frac{{\mathrm{\η}}_2^{\′}\·{\mathrm{\η}}_3^{\′}\·{\mathrm{\η}}_4^{\′}\·X(f_2,d_2)}{{\mathrm{\η}}_2\·{\mathrm{\η}}_3\·{\mathrm{\η}}_4\·X(f_2,d_1)}---\left(11\right)$

由(8)、(10)得，

Φ_out(f₂，d₂)-Φ_out(f₂，d₁)＝θ₂+θ₃(f₂)+θ₄+Φ(f₂，d₁)-θ′₂-θ′₃(f₂)-θ′₄-Φ(f₂，d₂)    (12)

5、由(5)式和(11)式，设参数A，令

$A=\frac{X(f_2,d_2)/X(f_2,d_1)}{X(f_1,d_2)/X(f_1,d_1)}=\frac{X(f_2,d_2)\·X(f_1,d_1)}{X(f_2,d_1)\·X(f_1,d_2)}=\frac{A{C}_{\mathrm{out}}(f_2,d_2)\·A{C}_{\mathrm{out}}(f_1,d_1)}{{\mathrm{AC}}_{\mathrm{out}}(f_2,d_1)\·{\mathrm{AC}}_{\mathrm{out}}(f_1,d_2)}---\left(13\right)$

由(6)式和(12)式，设参数Θ，令

Θ＝Φ(f₂，d₁)+Φ(f₁，d₂)-Φ(f₂，d₂)-Φ(f₁，d₁)

                                               (14)

＝[Φ_out(f₂，d₂)+Φ_out(f₁，d₁)-Φ_out(f₂，d₁)-Φ_out(f₁，d₂)]+[θ₃(f₂)+θ₃(f₁)-θ₃(f₁)-θ₃(f₂)]

从式(13)可以看出，A值仅与由受检体本身引起的幅度衰减X有关，实际中可以通过测量AC_out得到，不受系统固有因素的影响；式(14)中的Θ也仅与受检体本身引起的相位延迟Φ有关，但实际测量结果中仍包含系统固有因素θ₃和θ′₃的影响。在特定条件下，系统固有因素产生的相位延迟为一常量，利用标准体可测量得出其大小，设此常量参数Φ_cal，令

Φ_cal＝θ₃(f₂)+θ′₃(f₁)-θ₃(f₁)-θ′₃(f₂)         (15)

下面具体说明本发明的频域信息获取方法和时域参数重构方法。

1、利用频域蒙特卡洛模拟得到光学参数已知的标准体在光源调制频率为f₁、探测距离为d₁时的理论相位延迟

以及探测距离为d₂时的理论相位延迟

利用频域蒙特卡洛模拟得到光学参数已知的标准体在光源调制频率为f₂、探测距离为d₁时的理论相位延迟

以及探测距离为d₂时的理论相位延迟

(本发明中，采用MonteCarlo模拟的方法建立光学参数与实际中测量到的样品的光学特性之间的映射关系。蒙特卡洛计算模型的基本框架，主要包括光子随机步长的确定、光子吸收和散射过程的描述、以及待测物理量的统计描述。有关该方法的具体实现，可参见下文：MCML-Monte Carlo modelingof light transport in multi-layered tissues/Lihong Wang，Steven L.Jacques，Liqiong Zheng//Computer Methods and Programs in Biomedicine，Volume 47，Issue 2，July 1995，Pages131-146)。具体方法如下：

在两个探测频率和两个探测距离下，对一定范围(0.10cm^-1≤μ_a≤0.80cm^-1，30cm^-1≤μ′_s≤100cm^-1)等间隔的225组光学参数进行时域蒙特卡洛模拟。并利用快速公式结合朗伯比尔-定理进行时-频变换，见(16)、(17)式，分别得到各组光学参数对应的幅值AC和相位Φ。即建立4个225组光学参数下的模拟绝对频域信息数据库。

$u_{\mathrm{sc}}^{\%}\left({\mathrm{\ω}}_0\right)=\underset{k}{\mathrm{\Σ}}{W}_k\exp \left({\mathrm{i\ω}}_0{t}_k\right)---\left(16\right)$

W′_k＝W_kexp(μ_al_k)                                        (17)

式(19)为快速公式，式中

的实部和虚部分别是频域信息幅值AC和相位Φ，W_k表示第k步的光子包权重，ω₀为探测时的调制频率，t_k表示第k步的光子平均飞行时间；式(20)结合了朗伯-比尔定理，式中W′_k表示在μ_a≠0时的光子包权重，l_k表示第k步的光子平均距离。

其次，根据(13)、(14)式，分别将每一组光学参数在两个距离两个频率下的绝对幅值信息AC两次作比，绝对相位信息Φ两次做差，得到与测量量相对应的模拟相对频域信息A_th、Θ_th数据库。

2、令受检体为标准体，执行在前述的“频域参数获取的理论依据”里提到的探测步骤1-4，由式(14)、(15)计算得，

3、令受检体为未知组织体，执行在前述的“频域参数获取的理论依据”里提到的探测步骤1-4，由式(13)可得，

由式(14)、(16)得，

4、根据式(19)、(21)中A^组织体和Θ^组织体的大小，重构出未知组织体的光学参数，即吸收参数μ_a和漫反射参数μ′_s，过程如下(更为具体的过程可参见申请号为201010171262.6的发明“基于频域近红外光测量的光学参数重构方法”)：

第一，利用二元多项式分别对A_th和Θ_th构成的数据曲面进行拟合，建立数学模型，见(21)、(22)式，即μ_a和μ′_s到A_th、Θ_th的映射关系，以解决在重构的迭代中快速得到任意光学参数下的A_th和Θ_th。

$A_{\mathrm{th}}({\mathrm{\μ}}_a,{\mathrm{\μ}}_s^{\′},n)=(a_0+a_1{\mathrm{\μ}}_a+a_2{\mathrm{\μ}}_a^2+...+a_n{\mathrm{\μ}}_a^n)\×(b_0+b_1{\mathrm{\μ}}_s^{\′}+b_2{\mathrm{\μ}}_s^{\′2}+...+b_n{\mathrm{\μ}}_s^{\′n})---\left(21\right)$

${\mathrm{\Θ}}_{\mathrm{th}}({\mathrm{\μ}}_a,{\mathrm{\μ}}_s^{\′},n)=(c_0+c_1{\mathrm{\μ}}_a+c_2{\mathrm{\μ}}_a^2+...+c_n{\mathrm{\μ}}_a^n)\×(d_0+d_1{\mathrm{\μ}}_s^{\′}+d_2{\mathrm{\μ}}_s^{\′2}+...+d_n{\mathrm{\μ}}_s^{\′n})---\left(22\right)$

式中a，b，c，d各系数由最小二乘回归决定，n为二元多项式的阶数。

第二，以上二元多项式模型结合L-M非线性优化算法，即构成逆蒙特卡洛模拟的重构程序

最后，将测量得到的A^组织体和Θ^组织体代入根据上述方法编制的重构程序，选取适当的初值A₀和Θ₀，比较计算出新的迭代值，并利用二元多项式映射模型，快速得到更新后迭代值的A_th和Θ_th，再与A^组织体和Θ^组织体相比较，更新后继续循环，多次调用二元多项式模型，直至测量量A^组织体、Θ^组织体与模拟量A_th、Θ_th的平方误差达到最小，即可认为该迭代值为相应的待测光学参数，即吸收参数和得反射参数。

Claims (1)

1.一种近红外漫射光频域信息获取方法，包括下列步骤：

(1)首先建立一个双距离探测频域系统，在该系统中，光源为频率可调的近红外光源，用于为系统提供入射到被检组织的光源，光源发出的光分成两束，分别被导入到受检体的两个不同位置，被受检体在探测点处的漫反射光由探测器检测出，设两束光入射到受检体的两个不同位置距离探测点的距离分别为探测距离d₁和探测距离d₂，设系统的输入幅度为AC_in，初始相位为Φ_in，输出幅度为AC_out，输出相位为Φ_out；

(2)利用频域蒙特卡洛模拟得到光学参数已知的标准体在光源调制频率为f₁、探测距离为d₁时的理论相位

(f₁，d₁)以及探测距离为d₂时的理论相位

(f₁，d₂)；利用频域蒙特卡洛模拟得到光学参数已知的标准体在光源调制频率为f₂、探测距离为d₁时的理论相位

(f₂，d₁)以及探测距离为d₂时的理论相位(f₂，d₂)，并根据下列公式计算系统频域相位延迟参数Φ_cal：

(3)令受检体为标准体，分别探测光源调制频率为f₁、探测距离为d₁时标准体的输出相位

光源调制频率为f₁、探测距离为d₂时标准体的输出相位

探测光源调制频率为f₂、探测距离为d₁时标准体的输出相位

和探测光源调制频率为f₂、探测距离为d₂时的标准体输出相位

(4)令受检体为未知组织体，分别探测光源调制频率为f₁、探测距离为d₁时的输出幅度

(f₁，d₁)和输出相位

(f₁，d₁)，光源调制频率为f₁、探测距离为d₂时的输出幅度

(f₁，d₂)和输出相位

(f₁，d₂)，探测光源调制频率为f₂、探测距离为d₁时的输出幅度

(f₂，d₁)和输出相位

(f₂，d₁)以及探测光源调制频率为f₂、探测距离为d₂时的输出幅度

(f₂，d₂)和输出相位

(f₂，d₂)；

(5)根据下列公式计算未知组织体在频域中的幅值信息A^组织体：

(6)根据下列公式计算未知组织体在频域中的相位信息Θ^组织体：

Images