为克服现有技术的不足，鉴于生物组织对光的强散射作用及光子扩散过程的时空统一性， 其表面任意点的溢出光流均包含其内部荧光参数的空间分布信息，并反映在所测量光强的时 空变化模式上，由此，本发明的目的在于：提供一种面向小动物分子成像的时域荧光扩散层 析系统，获得高定量精度和空间分辨率的三维层析方法，以及稳健的双参数和多组分成像能 力。
本发明采用的技术方案是：面向小动物分子成像的时域荧光扩散层析系统，包括：
提供所需波长超短激光的两个皮秒半导体激光器及荧光探针；
用于选择前述两个皮秒半导体激光器出射的两路激光的一个2:1光开关；
用于将2:1光开关出射激光强度衰减到所需测量范围的可变衰减器；
用于将可变衰减器出射的所选波长光源切换成多个光源点的光开关；
用于将光源点出射激光投射到目标体的入射光纤；
用于接收来自目标体的反射或透射激光的接收光纤；
用于接收来自接收光纤的检测模块，检验模块进一步依次包括：滤除相应激发光的带通 滤波器、PMT光电倍增管检测器、TCSPC单光子计数模块；
用于由检验模块的输出产生测量平面或三维荧光发射率和寿命的图像的荧光分子层析时 域FMT图像重建模块，FMT图像重建模块进一步包括下列模块：
①对检测量Гx(ζd，ξs，t)和Гm(ζd，ξs，t)进行拉普拉斯变换，得到Γm(ζd，ξs，p)和Γx(ζd，ξs，p)的 模块；
②求出荧光波长检测量相对于激发波长检测量的Born ratio玻恩比： $R_{m/x}({\zeta }_d,{\xi }_s,p)=\frac{{\Gamma }_m({\zeta }_d,{\xi }_s,p)}{{\Gamma }_x({\zeta }_d,{\xi }_s,p)}$的模块；
③在Robin边界条件下，对激发光波长的拉普拉斯变换时域扩散方程求解，获得ξs处激励下r 和ζd处之光密度拉普拉斯变换值Φx(r，ξs，p)和Φx(ζd，ξs，p)的模块，即求解
$\left\{\begin{array}{c}[{▿·D}_x\left(r\right)▿-{\mu }_{\mathrm{ax}}\left(r\right)c-p]{\Phi }_x(r,{\xi }_s,p)=-\delta (r-{\xi }_s)\\ c{\Phi }_x(r,{\xi }_s,p)+2K{D}_x\left(r\right)n{·▿\Phi }_x(r,{\xi }_s,p){|}_{r\in \partial \Omega }\\ {\Phi }_x({\zeta }_d,{\xi }_s,p)={\Phi }_x(r,{\xi }_s,p){|}_{r={\zeta }_d}\end{array}=0\right.$
其中μax和Dx为激发波长下的吸收系数和扩散系数；Ω为成像域Ω之边界；n为边界检测点 之外法向；K为内反射控制因子；c为组织中的光速；
④在罗宾(Robin)边界条件下，获得荧光波长下的拉普拉斯变换时域扩散方程之格林函数 Gm(ζd，r，p)的模块；即求解
$\left\{\begin{array}{c}[{▿·D}_m\left(r^{\prime }\right)▿-{\mu }_{\mathrm{am}}\left(r^{\prime }\right)c-p]G_m(r^{\prime },r,p)=-\delta (r^{\prime }-r)\\ c{G}_m(r^{\prime },r,p)+2K{D}_m\left(r^{\prime }\right)n{·▿G}_m(r^{\prime },r,p){|}_{r^{\prime }\in \partial \Omega }\\ G_m({\zeta }_d,r,p)=(\mathrm{cK}/2)G_m(r^{\prime },r,p){|}_{r^{\prime }={\zeta }_d}\end{array}=0\right.$
其中r和r′为两个独立的空间矢量，μam和Dm为荧光波长下的吸收系数和扩散系数，Gm(r′，r，p)为 r处激励r′处出射的格林函数拉普拉斯变换值；
⑤由荧光扩散方程
$\left\{\begin{array}{c}[{▿·D}_m\left(r\right)▿-{\mu }_{\mathrm{am}}\left(r\right)c-p]{\Phi }_m(r,r_s,p)=-c{\Phi }_x(r,r_s,p)x(r,p)\\ x(r,p)={\mathrm{\eta \mu }}_{\mathrm{af}}\left(r\right)/[1+\mathrm{p\tau }\left(r\right)]\end{array}\right.$
得出荧光扩散成像积分方程Γm(ξd，ζs，p)＝∫ΩGm(ξd，r，p)Φx(r，ζs，p)x(r，p)dΩ，并进行空间 离散化获得成像矩阵方程Γm(p)＝W(p)x(p)的模块，其中x(p)＝[x1(p)，x2(p)，...，xN(p)]T； Γm(p)＝[Γm(ζ1，ξ1，p)，Гm(ζ2，ξ1，p)，...，Гm(ζD，ξS，p)]T；W(p)为SD×N维矩阵；N为离散后之体 元数；
⑥由任务②和⑤得Born比形式成像矩阵方程${\hat{\Gamma }}_m\left(p\right)=W\left(p\right)x\left(p\right)$之模块，其中中元素的 计算式为${\hat{\Gamma }}_m({\zeta }_d,{\xi }_s,p)=R_{m/x}({\zeta }_d,{\xi }_s,p)F^{\left(\Gamma \right)}({\zeta }_d,{\xi }_s,p),$F(Γ)(ζd，ξs，p)模型计算值；
⑦基于代数重建技术(Algebraic Reconstruction Technique，ART)的成像矩阵方程求解模块， 即
$\left\{\begin{array}{c}{x}_{k+1}\left(p\right)=x_k\left(p\right)+\lambda \frac{[{{\hat{\Gamma }}_m}^{(k+1)}\left(p\right)-W^{(k+1)}\left(p\right)x_k\left(p\right)]}{\left[W^{(k+1)}\left(p\right)\right]·{\left[W^{(k+1)}\left(p\right)\right]}^T}{\left[W^{(k+1)}\left(p\right)\right]}^T\\ k=\mathrm{0,1,2},...,(S\times D-1)\end{array}\right.$
其中为的第k个元素；W(k)(p)为W(p)的第k行；
⑦基于变换因子对p1和p2，以及上述ART方法的荧光产率ημaf(r)和时间寿命τ(r)成像计算模 块，即
$\left\{\begin{array}{c}{\mathrm{\eta \mu }}_{\mathrm{af}}\left(r\right)=(p_1-p_2)x(r,p_1)x(r,p_2)/[p_{1}x(r,p_1)-p_{2}x(r,p_2)]\\ \tau \left(r\right)=-[x(r,p_1)-x(r,p_2)]/[p_{1}x(r,p_1)-p_{2}x(r,p_2)]\end{array}\right.$
上述公式中r和r′为两个独立的空间矢量，ξs代表光源的位置，ζd为探测器的位置，μax和Dx为 激发波长下的吸收系数和扩散系数，Ω为成像域Ω之边界，n为边界检测点之外法向，K为 内反射控制因子，c为组织中的光速，μam和Dm为荧光波长下的吸收系数和扩散系数， Gm(r′，r，p)为r处激励r′处出射的格林函数拉普拉斯变换值，Γx(ζd，ξs，t)为光源ξs激励，ζd 点探测时激发光的检测光流量，Γx(ζd，ξs，p)为其拉普拉斯变换；Γm(ζd，ξs，t)为光源ξs激励， ζd点探测时出射光的检测光流量，Γm(ζd，ξs，p)为其拉普拉斯变换；Φx(r，ξs，p)为在外推边 界条件下光源ξs激励，r位置处的激发光的光密度；Φx(ζd，ξs，p)为在外推边界条件下光源ξs 激励，ζd位置处探测到的激发光的光密度；Gm(ζd，r，p)为r位置产生荧光时ζd位置处探测到 的出射光光流量的格林函数解；μav(r)，μsv′(r)，Kv(r，t)＝c/[3μsv′(r)]分别为吸收、退化、 散射和扩散系数，其中v为x时表示激发光，v为m时表示出射光；ημaf(r)表示荧光产率；η 为量子效率；μaf(r)为荧光吸收系数；τ(r)为荧光寿命，δ(r，rs)为rs处的弥向点光源， Φx(r，rs，p)为rs处光源激励，r出激发光的光密度，Φm(r，rs，p)为rs处光源激励，r出激发光 的光密度。
所述的入射光纤2和同轴混合光纤3分布在目标体的同侧，目标体被压缩成平板，所 述切换成多个光源点的光开关为切换成十六个，形成十六个光源光纤入射点，在目标体的 上表面分四行四列布置十六个光源光纤入射点，其中中间四个点的光纤由同轴混合光纤3 组成，每个同轴混合光纤的接收光纤束5分别连接四个带通滤波器中的一个，每个带通滤波 器又通过四个PMT检测器中的一个连接到四通道TCSPC模块的一个通道，所述的求出外 推边界条件下时域扩散方程的拉普拉斯变换后的解Φx(r，ξs，p)、Φx(ζd，ξs，p)和Gm(ζd，r，p) 的模块是利用其外推边界条件下的时域扩散方程的拉普拉斯变换，可分别求出其半无限条 件下的解析解。
所述的入射光纤2和接收光纤6分布在目标体的两侧，所述切换成多个光源点的光开 关为切换成十六个，形成十六个光源光纤入射点，目标体的入射表面以四行四列布置十六 个光源光纤入射点，其探测表面以两行两列布置四个接收光纤，每个接收光纤各自通过四 个PMT检测器中的一个连接到四通道TCSPC模块的一个通道，所述的求出外推边界条件 下时域扩散方程的拉普拉斯变换后的解Φx(r，ξs，p)、Φx(ζd，ξs，p)和Gm(ζd，r，p)的模块是利 用外推边界条件下的时域扩散方程的拉普拉斯变换，可分别求出其半无限条件下的解析解。
所述的光纤由同轴混合光纤3组成，六十四个同轴混合光纤3围绕目标体呈四层分布， 每层沿圆周方向均匀分布十六个同轴混合光纤3作为检测点，所述切换成多个光源点的光 开关为切换成十六个，形成十六个光源光纤入射点，每次连接一层，四个4:1光开关分四 次、每次选择目标体上十六个检测点中的四个检测点对应连接到四个带通滤波器，每个带 通滤波器分别对应连接四个PMT检测器中的一个，并通过PMT检测器传送至四通道TCSPC 模块，所述的求出外推边界条件下时域扩散方程的拉普拉斯变换后的解Φx(r，ξs，p)、 Φx(ζd，ξs，p)和Gm(ζd，r，p)的模块是利用有限元方法求出其时域扩散方程的拉普拉斯变换 的数值解。
所述两个皮秒半导体激光器的波长依据选用的荧光探针有两种组合方式：选用波长为 785nm和830nm，配合的荧光探针为ICG；选用波长为670nm和710nm，配合的荧光探针 为CY5.5。
本发明基于时间分辨测量，提供了一种面向小动物分子成像的时域荧光扩散层析方法与 系统，不仅可以可有效地涵盖荧光发射的时空分布信息，并从本质上提供了两参数(荧光产 率和寿命)、多组分荧光图像同时重建能力，进而实现活体小动物FMT成像功能和质量的实 质性提高。采用这种面向小动物分子成像的时域荧光扩散层析方法与系统，使活体生物内部 的特异分子生化反应及其变化过程的实时定量观测成为现实，对阐释生命活动规律，揭示疾 病发生机理，建立疾病预警机制，提高诊疗水平和创制新药物都具有重大意义。
附图说明
图1多通道时间分辩FMT实验系统的基本构成框图。
图2平板压缩反射式检测模型示意图。
图3同轴混合光纤束探头的结构示意图。
图4平板压缩透射式检测模型示意图。
图5全三维圆柱压缩透射式的检测模块结构框图。
图6平板压缩反射式(或透射式)的检测模块结构框图。
图7全三维圆柱压缩透射式检测模型示意图。
图8荧光分子层析图像重建方法流程图。
图中，1为肿瘤目标，2为入射光纤，3为同轴混合光纤，4为同轴混合光纤的入射光纤， 5为同轴混合光纤的接收光纤束，6为接收光纤。

本发明包括以下两部分：多通道时域FMT测量系统；基于广义脉冲谱技术(Generalized Pulse Spectrum Technique，GPST)的荧光分子层析图像重建方法。
A.时域荧光分子层析测量技术与系统：
时域FMT技术研究旨在发展高灵敏时间分辨检测技术，通过多“角度”激发测量表面荧光 瞬态“投影”和基于精确光子迁移模型的图像反演算法，重建组织体内部特定分子靶标荧光参 数的空间分布，并通过时间分辨信息弥补空间采样数量的有限性以有效提高成像质量和检测 灵敏度，实现FMT的应用要求。
皮秒或飞秒超短激光脉冲通过高散射组织体后的光子时间分布称为时间点扩展函数 (Temporal Point Spread Function，TPSF)，其幅度和形状蕴涵了组织体内部光学参数分布的 信息。对于FMT，上述物理现象涉及两个耦合的光子迁移过程，其中激发光在组织体内扩散 传播时与荧光探针相互作用(吸收)而形成二次荧光发射，其时空特性不仅取决于荧光探针的 发射率和寿命，而且与激发光的扩散强度正比相关，因此，建立组织体内荧光产生和迁移过 程的数学模型并测量表面荧光发射的TPSF即可还原其内部的荧光参数分布。
随着时间相关单光子计数(TCSPC)和超短激光脉冲技术的进步，时间分辨测量的成本 和复杂性大为降低，基于数字式TCSPC的多通道时间分辨系统成为DOT等散射光成像研究的 主要测量技术。TCSPC的基本原理是通过测量单光子到达时间与固定参考时间之差，建立与 TPSF成正比的光子飞行时间统计直方图。该技术的时间分辨性能主要取决于光电倍增管 (PMT)检测器的渡越时间弥散(Transient Time Spread，TTS)特性。目前电极式PMT的TTS 最小可达约150ps，而微通道板PMT则可实现小于50ps的TTS。TCSPC结合了光子计数和超 快电子技术的综合优势，具有灵敏度高、动态范围大以及时间分辨率合理等一系列突出优点， 因而非常适合时域FMT中具有纳秒持续时间量级的超微弱瞬态荧光信号的检测。在工作原理 上，为防止光子“堆积”效应，TCSPC的最高计数率一般不超过光源重复频率的1/10，因此基 于TCSPC原理的时间分辨测量系统在设计上必须求得计数率(测量时间)和量程间的平衡。
鉴于大多数红外荧光探针的衰减时间为0.5-4ns，而在基于小动物模式FMT和一般荧光层 析技术研究和验证中，源-探测器最大间距约为3-5厘米，因此可得激发和荧光输出光流时间 点扩展函数的时间范围大约为2-6ns(计算基于目前公认的软组织光学参数和普遍应用的荧光 探针ICG)。故本研究中所需要的测量时间范围不小于10ns，同时系统应具有较高采样分辨 率和较短的测量时间，以便使FMT的应用能够有效反映生化过程的时间变化规律。
根据以上基本考虑和TCSPC技术的最新技术发展，本发明构建的多通道高速TCSPC时间 分辨FMT测量系统的基本构成如图1所示，由皮秒半导体激光器、2:1光开关、可变衰减器、 1:16光开关、目标体和检测模块组成。
●皮秒半导体激光器提供所需波长的超短激光光源；
●2:1光开关用于选择皮秒半导体激光器出射的两路光源；
●可变衰减器用于将光源强度衰减到所需测量范围；
●1:16光开关用于将所选波长光源切换成16个光源点；
●目标体依据测量方式不同，有两种处理方式：一种是对目标体进行平板压缩，该方式 应用于平板压缩反射式(见图2)或是平板压缩透射式(见图4)测量方式；另一种是 全三维圆柱压缩透射式(见图5)，该方式将目标体放置于圆柱形成像腔中，并用组织 模拟液对成像腔的空隙部分进行填充，该方式适用于三维透射测量。在实施例中目标 体为荷瘤小动物(实验鼠)。
●检测模块依据检测方式不同，分为四通道检测和切换式十六通道检测两种方式。其中 前者主要用于平板压缩反射式或平板压缩透射式测量方式，该检测方式的检测模块结 构框图如图6所示，对于每次源入射，来自目标体的四个检测点分别经过带通或长通 滤光器，滤除相应激发光，然后分别由PMT检测器进行检测，并传送至四通道TCSPC 模块；后者用于三维透射测量，该检测方式的检测模块结构框图如图7所示，对于每 次源入射，四个4:1光开关分四次、每次选择被测目标体成像面上十六个检测点中的 四个检测点作四通道并行TCSPC检测，通过分别调节四个光路的带通或长通滤光器滤 除相应激发光，然后分别由PMT检测器进行荧光信号检测，并传送至四通道TCSPC 模块，实现切换式十六通道的检测；
本发明针对平板压缩反射式检测模式，设计的检测模型如图2所示，入射光纤2和同轴混 合光纤3分布在目标体的同侧，在目标体的上表面分四行四列布置了十六个光源光纤入射点， 其中中间四个点的光纤由同轴混合光纤3组成，该探头的结构如图3所示，中间为同轴混合光 纤的入射光纤4，周围为同轴混合光纤的接受光纤束5。该同轴混合光纤3可以精确地获得测量 信号并有效减小位于组织体上的探头数目。
本发明针对平板压缩透射式检测模式，设计的检测模型如图4所示，入射光纤2和接收光 纤6分布在目标体的两侧，目标体的入射表面以四行四列布置了十六个光源光纤入射点，其 探测表面以两行两列布置了四个探测光纤。该光源及探测方式的组合既能很好的获得目标体 的内部信息，也能有效的减少成本，具有很好的性能价格比。
本发明针对全三维圆柱压缩透射式检测模式，设计的检测模型如图5所示，光纤探头仍 采用图5所示同轴混合光纤3，六十四个光纤探测点围绕目标体呈四层分布，每层沿圆周方 向均匀分布十六个光纤探头点。每次源入射，与光源同层的十六个检测点作并行TCSPC检测， 检测结束后切换源位置重复上述测量，直至本层十六个源位置“扫描”完毕并转至下一层，重 复上述过程，直到四层全部扫描完毕，最后可利用时域FMT图像重建算法产生测量平面或全 三维荧光产率和寿命的图像。该光源及探测方式的组合既能很好的获得目标体的内部信息， 也能有效的减少成本，具有很好的性能价格比。
B.基于GPST的荧光分子层析图像重建理论和方法
图像重建算法是FMT技术最重要的研究课题之一，它直接关系到成像系统方案的实现和 所能达到的性能指标。FMT反演方法研究必须以有效实现实际应用所要求的高分辨率、高精 度、通用性、实时性和稳健性为宗旨，而时域测量模式为该目标的实现提供了有力的手段和 技术可行性。
在FMT的实现中，由于组织体的强散射效应和解剖结构不均匀性，激发光和荧光均呈现 复杂的传播模式，其数学模型可近似地由一对耦合光子扩散方程描述。其中激发光在漫射传 播中作为荧光激发源，因此组织体各点荧光发射强度不但正比于入射激发光强度，而且取决 于激发波长下的光学参数，其分布图像强烈地影响激发光子密度的分布，进而相当程度地影 响各检测点荧光发射强度及表面时域测量的瞬态特性。实现具有高特异性和灵敏度的时域 FMT的完善方案应包括三个重建过程：即首先由一般DOT过程分别实现激发光和发射光波长 下的光学参数的重建，及在此基础通过组织体内部受激荧光发射的表面测量模式重建组织内 部的荧光参数的空间分布。
遵循此原则，并基于牛顿-拉夫逊框架的时域DOT图像重建算法，本发明提出了一个频率 对条件下的时域GPST-FMT单组分、两参数图像同时重建算法，其流程如图8所示。
假设ξs(s＝1，2，...，S)代表光源的位置，ζd(d＝1，2，...，D)为探测器的位置，r为荧光产生 的位置，μav(r)，μsv′(r)，Kv(r，t)＝c/[3μsv′(r)]分别为吸收、退化、散射和扩散系数，其 中v为x时表示激发光，v为m时表示出射光；c为光在组织中的速度，ημaf(r)表示荧光产率； η为量子效率；μaf(r)为荧光吸收系数；τ(r)为荧光寿命，δ(r，rs)为rs处的弥向点光源，则 时域荧光耦合扩散方程进行拉普拉斯变换得
$\left\{\begin{array}{c}[{▿·K}_x\left(r\right)▿-{\mu }_{\mathrm{ax}}\left(r\right)c-p]{\Phi }_x(r,r_s,p)=-\delta (r,r_s)\\ [{▿·K}_m\left(r\right)▿-{\mu }_{\mathrm{am}}\left(r\right)c-p]{\Phi }_m(r,r_s,p)=-c{\Phi }_x(r,r_s,p)x(r,p)\\ x(r,p)={\mathrm{\eta \mu }}_{\mathrm{af}}\left(r\right)/[1+\mathrm{p\tau }\left(r\right)]\end{array}\right.---\left(1\right)$
其中Φx(r，rs，p)为rs处光源激励，r出激发光的光密度，Φm(r，rs，p)为rs处光源激励，r出激 发光的光密度。
针对平板压缩反射式或是透射式均匀检测模型，利用其外推边界条件下的时域扩散方程 的拉普拉斯变换，可分别求出其半无限条件下的解析解Φx(r，ξs，p)、Φx(ζd，ξs，p)和 Gm(ζd，r，p)；而对于三维圆柱压缩透射式检测模型，可利用有限元方法求出其时域扩散方程 的拉普拉斯变换的数值解Φx(r，ξs，p)、Φx(ζd，ξs，p)和Gm(ζd，r，p)。 假设Γx(ζd，ξs，t)为光源ξs点激励、ζd点探测时激发光的检测光流量，Γx(ζd，ξs，p)为其拉普 拉斯变换；Γm(ζd，ξs，t)为光源在ξs点激励、ζd点探测时出射光的检测光流量，Γm(ζd，ξs，p) 为其拉普拉斯变换。根据Born Ratio，有
$R_{m/x}({\zeta }_d,{\xi }_s,p)=\frac{{\Gamma }_m({\zeta }_d,{\xi }_s,p)}{{\Gamma }_x({\zeta }_d,{\xi }_s,p)}---\left(2\right)$
由方程(1)之第二式可推出${\Gamma }_m({\xi }_d,{\zeta }_s,p)=\underset{\Omega }{\int }c{G}_m({\xi }_d,r,p){\Phi }_x(r,{\zeta }_s,p)x\left(r\right)\mathrm{d\Omega },$对其进行离散 化，可得
          Γm(p)＝W(p)x(p)            (3)
其中x(p)＝[x1(p)，x2(p)，...，xN(p)]T；Γm(p)＝[Гm(ζ1，ξ1，p)，Гm(ζ2，ξ1，p)，...，Гm(ζD，ξS，p)]T；W(p) 为SD×N维矩阵；N为离散后之体元数。W(p)的元素计算如下：
          W(ζd，ξs，n，p)＝cGm(ζd，rn，p)Φx(rn，ξs，p)Vn     (4)
其中rn为第n个体元的中心位矢，Vn为第n个体元的体积。
方程(3)玻恩比形式为
${\hat{\Gamma }}_m\left(p\right)=W\left(p\right)x\left(p\right)---\left(5\right)$
其中中元素的计算式为Гm(ζd，ξs，p)＝Гm(ζd，ξs，p)·F(Γ)(ζd，ξs，p)，F(Γ)(ζd，ξs，p)模型计 算量。
利用代数重建技术(ART)求解上式中的x(p)，有
$\left\{\begin{array}{c}{x}_{k+1}\left(p\right)=x_k\left(p\right)+\lambda \frac{[{{\hat{\Gamma }}_m}^{(k+1)}\left(p\right)-W^{(k+1)}\left(p\right)x_k\left(p\right)]}{\left[W^{(k+1)}\left(p\right)\right]·{\left[W^{(k+1)}\left(p\right)\right]}^T}{\left[W^{(k+1)}\left(p\right)\right]}^T\\ k=\mathrm{0,1,2},...,(S\times D-1)\end{array}\right.---\left(6\right)$
其中为的第k个元素；W(k)(p)为W(p)的第k行；λ为迭代松弛因子。
将拉普拉斯变换因子p分别取p1和p2，利用ART技术求解式(5)可分别得到x(p1)和x(p2)， 依据方程(1)之第三式，可进一部求出荧光产率ημaf(r)和时间寿命τ(r)
$\left\{\begin{array}{c}{\mathrm{\eta \mu }}_{\mathrm{af}}\left(r\right)=(p_1-p_2)x(r,p_1)x(r,p_2)/[p_{1}x(r,p_1)-p_{2}x(r,p_2)]\\ \tau \left(r\right)=-[x(r,p_1)-x(r,p_2)]/[p_{1}x(r,p_1)-p_{2}x(r,p_2)]\end{array}\right.---\left(7\right)$
上述时域荧光分子层析图像重建算法以光子扩散理论为依据，理论严密，功能完备，其 主要特点是基于特征数据、计算速度快、数据使用灵活以及稳健性强，因而特别适合荧光扩 散层析中多参数三维重建问题。
下面结合附图和实施例对本发明做进一步的详细说明。
图1所示，多通道时间分辩FMT实验系统主要由皮秒半导体激光器、2:1光开关、可变 衰减器、1:16光开关、目标体和检测模块组成。其中皮秒半导体激光器的脉宽小于50ps、重 复频率为50MHz、平均功率小于1mW，依照所选荧光探针的不同可选用两组波长①皮秒半 导体激光器分别选用波长为785nm和830nm，配合的荧光探针为CY5.5②皮秒半导体激光器 分别选用波长为670nm和710nm，配合的荧光探针为ICG。该系统的基本工作过程设定如下： 入射光纤2、同轴混合光纤3和接收光纤6根据应用需要交叉布置在成像平面或呈三维配置， 2:1光开关用于选择不同波长光源，可变衰减器将所选波长光源的强度衰减到所需测量范围， 1:16光开关将选定的光源由入射光纤或是同轴混合光纤的入射光纤4依次导入目标体的十六 个源位置。目标体依据检测方式不同，可分为平板压缩反射式检测模型、平板压缩透射式检 测模型和全三维圆柱压缩透射式检测模型，其模型示意图分别如图2、图4和图5所示。检 测模块用于检测组织表面的光流量，其具体检测方式依据检测模型的不同分为平板压缩反射 式(或透射式)的检测模块和全三维圆柱压缩透射式的检测模块，其具体结构框图见图6和 图7。
图2给出了平面压缩反射式检测模型示意图。目标体被压缩形成平面，入射光纤2与同 轴混合光纤3放置在目标体的平面一侧，实现反射式测量。十六个光源光纤以四行四列方式 布置，其中中心位置处由四个同轴混合光纤3组成，其具体结构参见图3。对每次源入射， 四个检测点作并行TCSPC检测，检测结束后切换源位置重复上述测量，直至十六个源位置“扫 描”完毕。
图3给出同轴混合光纤探测头的结构实施图。探测头中心位置为同轴混合光纤的入射光纤 4，周围为同轴混合光纤的接收光纤束5。该同轴混合光纤3探头可以精确地获得测量信号并有 效减小位于组织体上的探头数目。
图4给出了平板压缩透射式检测模型示意图，其中图(a)为十六个入射光纤排列示意图， 图(b)为四个接收光纤示意图。入射光纤2和接收光纤6分布在目标体的两侧，目标体的入射 表面以四行四列布置了十六个入射光纤2，其探测表面以两行两列布置了四个接收光纤6。对 每次源入射，四个检测点作并行TCSPC检测，检测结束后切换源位置重复上述测量，直至十 六个源位置“扫描”完毕。该光源及探测方式的组合既能很好的获得目标体的内部信息，也能 有效的减少成本，具有很好的性能价格比。
图5给出了全三维圆柱形压缩透射式检测模型示意图。光纤探头仍采用图5所示同轴混 合光纤3，六十四个光纤探测点围绕目标体呈四层分布，每层沿圆周方向均匀分布十六个同 轴混合光纤3。对每次源入射，与光源同层的十六个检测点作并行TCSPC检测，检测结束后 切换源位置重复上述测量，直至本层十六个源位置“扫描”完毕并转至下一层，重复上述过程， 直到四层全部扫描完毕，最后可利用时域FMT图像重建算法产生测量平面或三维荧光发射率 和寿命的图像。该光源及探测方式的组合既能很好的获得目标体的内部信息，也能有效的减 少成本，具有很好的性能价格比。
图6给出了平板压缩反射式(或透射式)的检测模块结构框图。对每次源入射，四个检 测点作并行TCSPC检测，分别调节四个光路的带通滤波器滤除相应激发光，保证PMT检测 器(PMT检测器内置了其工作所需高压电源和宽带放大器等，仪器响应小于200ps)工作在 可靠的单光子计数状态，四通道TCSPC模块拥有四个并行信号通道，其单通道的最高计数率 为8MHz、最小电子分辨为8ps。检测结束后切换源位置重复上述测量，直至十六个源位置“扫 描”完毕，最后可利用时域FMT图像重建算法产生测量平面或三维荧光发射率和寿命的图像。 上述十六源并行四通道检测方式兼顾测量时间和系统成本两方面的要求，具有很好的性能价 格比。
图7给出了全三维圆柱压缩透射式的检测模块结构框图。对每次源入射，四个4:1光开关 分四次、每次选择十六个检测点中的四个检测点作四通道并行TCSPC检测，通过分别调节四 个光路的带通滤波器滤除相应激发光，保证PMT检测器(PMT检测器内置了其工作所需高 压电源和宽带放大器等，仪器响应小于200ps)工作在可靠的单光子计数状态，之后切换源位 置重复上述测量，直至本层十六个源位置“扫描”完毕并转至下一层，重复上述过程，直到四 层六十四个探测点全部扫描完毕。四通道TCSPC模块拥有四个并行信号通道，其单通道的最 高计数率为8MHz、最小电子分辨为8ps，检测结束后切换源位置重复上述测量，最后可利 用时域FMT图像重建算法产生测量平面或三维荧光发射率和寿命的图像。上述六十四个源和 探测器分四层分布，采用并行四通道切换方式实现每层十六个检测点的检测方式既可兼顾测 量时间和系统成本两方面的要求，又可有效减小通道不一致性造成的系统标定难度和成像误 差，具有很好的性能价格比。
图8给出了荧光分子层析图像重建算法的流程图。该流程可分以下几个步骤：
①对检测量Γx(ζd，ξs，t)和Γm(ζd，ξs，t)进行拉普拉斯变换，得到Γm(ζd，ξs，p)和Γx(ζd，ξs，p)；
②求出Born ratio，$R_{m/x}({\zeta }_d,{\xi }_s,p)=\frac{{\Gamma }_m({\zeta }_d,{\xi }_s,p)}{{\Gamma }_x({\zeta }_d,{\xi }_s,p)};$
③对荧光耦合扩散方程进行拉普拉斯变换，求出外推边界条件下的Φx(r，ξs，p)、Φx(ζd，ξs，p) 和Gm(ζd，r，p)(针对平面压缩反射式或平板透射式检测模型，利用其外推边界条件下的时域 扩散方程的拉普拉斯变换，可分别求出其半无限条件下的解析解Φx(r，ξs，p)、Φx(ζd，ξs，p)和 Gm(ζd，r，p)；而对于三维圆柱压缩透射式检测模型，可利用有限元方法求出其时域扩散方程 的拉普拉斯变换的数值解Φx(r，ξs，p)、Φx(ζd，ξs，p)和Gm(ζd，r，p))。
④由荧光耦合扩散方程
$\left\{\begin{array}{c}[{▿·K}_x\left(r\right)▿-{\mu }_{\mathrm{ax}}\left(r\right)c-p]{\Phi }_x(r,r_s,p)=-\delta (r,r_s)\\ [{▿·K}_m\left(r\right)▿-{\mu }_{\mathrm{am}}\left(r\right)c-p]{\Phi }_m(r,r_s,p)=-c{\Phi }_x(r,r_s,p)x(r,p)\\ x(r,p)={\mathrm{\eta \mu }}_{\mathrm{af}}\left(r\right)/[1+\mathrm{p\tau }\left(r\right)]\end{array}\right.$
可推出积分方程${\Gamma }_m({\xi }_d,{\zeta }_s,p)=\underset{\Omega }{\int }c{G}_m({\xi }_d,r,p){\Phi }_x(r,{\zeta }_s,p)x\left(r\right)\mathrm{d\Omega },$对其进行离散化，可得成像 代数方程：Γm(p)＝W(p)x(p)
⑤由步骤②和④得波恩比形式的成像代数方程：${\hat{\Gamma }}_m\left(p\right)=W\left(p\right)x\left(p\right),$其中中元素的计 算式为Γm(ζd，ξs，p)＝Γm(ζd，ξs，p)·F(Γ)(ζd，ξs，p)，F(Γ)(ζd，ξs，p)模型计算量。
⑥用代数重建技术(ART)求解上式，得x(p)
$\left\{\begin{array}{c}{x}_{k+1}\left(p\right)=x_k\left(p\right)+\lambda \frac{[{{\hat{\Gamma }}_m}^{(k+1)}\left(p\right)-W^{(k+1)}\left(p\right)x_k\left(p\right)]}{\left[W^{(k+1)}\left(p\right)\right]·{\left[W^{(k+1)}\left(p\right)\right]}^T}{\left[W^{(k+1)}\left(p\right)\right]}^T\\ k=\mathrm{0,1,2},...,(S\times D-1)\end{array}\right.$
⑦取变换因子对p1和p2，利用ART技术求解式(5)分别的x(p1)和x(p2)，依据方程(1)，求出荧 光产率ημaf(r)和时间寿命τ(r)
$\left\{\begin{array}{c}{\mathrm{\eta \mu }}_{\mathrm{af}}\left(r\right)=(p_1-p_2)x(r,p_1)x(r,p_2)/[p_{1}x(r,p_1)-p_{2}x(r,p_2)]\\ \tau \left(r\right)=-[x(r,p_1)-x(r,p_2)]/[p_{1}x(r,p_1)-p_{2}x(r,p_2)]\end{array}\right.$