，/f)是相互独立的高斯白噪声，0(r，k)表示与药代动力学渗 透率相关的中间矢量：
J表示由第i间室向第j间 室的渗透率兩康示由第j间室向第i间室的渗透率，K康示由第i间室通过代谢系统排出体 外的排泄率，i = l，2,…，Nc;j = l，2,…，Ne;
[0029] 步骤=、药代动力学参数直接成像状态方程：
[0030] 结合全域内所有有限元剖分节点rn处的间室巧光剂浓度和间接药代动力学参数， n=l，2,…，N，定义待估计状态矢量X=[cT(；ri) ,0T(ri)，…，cT(rN)，0T(rN)]T;引入巧光层析 成像方程，得到直接成像状态方程，
[0031]
C終
[0032] 式(6)中，A = diag[K(ri)，I，K(r2)，I，...，K(rN)，I]，K(rn)由式(4)计算，I 为单位矩 阵，A化）=巧In(IO)W化化为测量矩阵，Tl和e分别是巧光剂的量子效率和消光系数，E是的于 (W,+Wt=)iV列的矩阵将X转换到总浓度Cr(二Ci+C;十…+ L\)，Ci=[Ci(ri)，…，Ci(rN)]T， W(k)是Nd行N列的巧光层析成像方程的权重矩阵，S化)和丫化)是独立的高斯白噪声，其协 方差矩阵为Q和R;
[0033] 步骤四、自适应扩展卡尔曼滤波是通过采用遗忘因子弥补对系统初始状态先验知 识缺失所带来的影响，基于自适应扩展卡尔曼滤波获取药代动力学参数图像：
[0034] 定义X(/c|/c-i)与P化Ik-I)分别为根据k-1步估计值预测的第k步的状态矢量和误 差协方差矩阵，x(/c…U与p(k-l)分别为k-1步的状态矢量和误差协方差矩阵的估计值， X(Af)和P化)分别为k步的状态矢量和误差协方差矩阵的估计值；
[0035] 具体步骤如下：
[0036] 预测步骤:根据上一步估计的状态矢量义(/C-1)与误差协方差矩阵P化-1)，预测当 前步的状态矢量义a jA，-1)与误差协方差矩阵P化I k-1)，
[0037]
(8)
[0038] 式(8)中，Mk)为遗忘因子，J化-1)为第k-1步的系统化cobin矩阵；
[0039] Kalman增益步骤:根据上述预测值计算当前步Kalman增益G化），并根据当前步测 量伯恩比Inb化)计算新息(1化），
[0040]
(12)
[0041] 更新步骤:根据当前步预测值与新息值更新当前步的状态矢量X作），误差协方差 矩阵P化)W及系统的化CObin矩阵J化），
[0042]
(巧）
[0043] 重复上述预测-Kalman增益-更新过程，直至得到状态矢量X的最终估计值义(/〇， 利用X=[cT(rl)，0T(rl)，…，cT(rN)，0T(rN)]T得到间室浓度和中间矢量0(r，k)，利用式(4)得 到渗透率Kij、渗透率Kji和排泄率Ki。
[0044] 与现有技术相比，本发明的有益效果是：
[0045] (1)本发明方法的直接重建方法的每一步估计不需要完整=维数据，而是采用螺 旋离散扫描方式获得即时并行测量数据，降低了对数据采集设备的高时间分辨率的要求， 实现成本较低的测量模式；
[0046] (2)本发明方法建立的多间室动力学模型可W更为准确地描述组织体，其参数具 有相应的实际生理意义，对癌症诊断具有较大潜力；
[0047] (3)本发明方法采用的直接重建状态方程去除了间接法的中间浓度重建过程，避 免了中间传统巧光层析成像中由于其逆问题病态性带来的误差。
[0048] (4)本发明方法采用自适应扩展卡尔曼滤波，通过采用遗忘因子弥补对系统初始 状态先验知识缺失所带来的影响，相比于传统扩展卡尔曼滤波实现更佳的重建结果。
【附图说明】
[0049] 图1是本发明所采用的多路复用探测平面示意图；
[0050] 图2是本发明中螺旋离散扫描数据组成示意图；
[0051 ]图3是本发明中螺旋扫描成像腔及光源位置示意图；
[0052] 图4是本发明实施例二室模型示意图。
【具体实施方式】
[0053] 下面结合附图和具体实施例对本发明技术方案作进一步详细描述，所描述的具体 实施例仅对本发明进行解释说明，并不用W限制本发明。
[0054] 本发明提出的一种螺旋离散扫描式巧光剂药代动力学参数直接成像方法，其中， 静脉注射巧光剂的组织体放置于一圆柱的成像腔中，包括W下步骤：
[0055] 步骤一、根据=维螺旋离散扫描模式获取巧光剂扫描数据：
[0056] 步骤1-1、光源经过准直后入射到成像腔表面，在成像腔的外周表面且与光源同水 平面的高度上（图1示出了假设光源为0%位于光源同一水平面(XY平面）101.25°到258.75° 的位置)均布D个探测位置妒（= 1，2,…，凸)，Nd个探测器W多路复用的方式测量，Nd卽， 采样间隔为AT;采用如上所述的配置，L次并行测量获得该光源下的所有探测位置的投影 数据，L = D/Nd;然后，成像腔旋转一个角度为A Angle,同时上升一个高度为A Z，W此类推， 直至完成全=维的扫描测量，在该扫描测量过程形成了S个光源位置rs(s = l，2,…，S);图3 所展示的成像腔半径为15mm，高为40mm，扫描从Zl = 36mm到Z2 = 5mm，Z轴方向步距A Z = 1mm，旋转角度AAngle = 22.5°，光源位置数8 = 32。
[0057] 步骤1-2重复P次步骤1-1所述的全=维的扫描测量，获取所需的扫描数据M(如何 获取该扫描数据M属于本技术领域的惯用技术手段，在此不再寶述），P次扫描测量所覆盖的 巧光剂的代谢时间为P X S X L X A T，如图2所示，
[0化引
（1)
[0059] 式（I)中，r (口，3) = [1。6(口，3，1)，1加（口，3,2)，...，1。6(口，3，1^)]了，其中，口=1, 2,…，P，s = l，2,…，S，l = l，2,…，L，p为S维扫描数据的序号，S为光源位置序号，1表示采 样序号，r (p，s)表示第P组=维扫描数据中的第S组投影数据、且由边界测量光流量的伯恩 比组成，
[0060]
C2)
[0061 ]式(2)中，Ix(rd，rs，t)和Im(rd，rs，t)分别是激发光和巧光的边界测量光流量，设Inb (P，S，1)为Inb化），即t = k A T时的并行测量数据，k=l，2,…，K，K = PXSXL
[0062] 步骤二、建立多间室巧光剂动力学模型：
[0063] 将成像腔内的组织体概念性的划分为N。个间室，利用多间室巧光剂动力学模型描 述不同间室之间的相互渗透W及通过代谢系统排出组织体外的过程，该多间室巧光剂动力 学模型如下：
[0064]
巧）
[0065] 式（3)中，巧|-,^^) = [('1(1*，/0.(：':(1',/:)，-"，〔\(|-.的'，。(1'，10表示第1个间室在时刻1^ A T时的间室浓度，i = 1，2，…，N。;K( 0 (r，k))是与药代动力学渗透率相关的系统矩阵，
[0066]
(4)
[0067] 式(4)中，《(r，k)和?(r，W是相互独立的高斯白噪声，9(r，k)表示与药代动力学渗 透率相关的中间矢量，
Ku表示由第i间室向第j间 室的渗透率，扣表示由第j间室向第i间室的渗透率，K康示由第i间室通过代谢系统排出体 外的排泄率，i=l，2,…，Nc;j=l，2,…，N。;常用二室模型将组织体分为血管部分和血管外部 分，并认为只有血管部分可通过身体循环系统把巧光剂排出体外，如图4所示，其系统矩阵 为
[006引
（5)
[0069] 其中，Kpe(r)与Kep(r)分别表示由血管部分渗透到血管外部分和由血管外部分渗 透回血管部分的渗透率，Kp为由血管部分直接排出体外的排泄率；
[0070] 步骤=、药代动力学参数直接成像状态方程：
[0071] 结合全域内所有有限元剖分节点rn处的间室巧光剂浓度和间接药代动力学参数， n=l，2,…，N，定义待估计状态矢量X=[cT(ri)，0T(ri)，…，cT(rN)，0T(rN)]T;引入巧光层析 DFT成像方程，得到直接成像状态方程，
[0072]
(約
[0073] 式(6)中，A = diag[K(;ri)，I，K(r2),K(TN)，1] ,K(Tn)由式(4)计算，I 为单位矩 阵，A化）=巧In(IO)W化化为测量矩阵，n和e分别是巧光剂的量子效率和消光系数，E是N行 仇+]〇~列的矩阵将蛾换到总浓度(：了(=<：1+<：2+'" +仁,^)，。=[(：加1)^''，。柄）]了， W(k)