CN112037297B - 一种基于网格重组的混合光传输建模方法 - Google Patents
一种基于网格重组的混合光传输建模方法 Download PDFInfo
- Publication number
- CN112037297B CN112037297B CN202010742538.5A CN202010742538A CN112037297B CN 112037297 B CN112037297 B CN 112037297B CN 202010742538 A CN202010742538 A CN 202010742538A CN 112037297 B CN112037297 B CN 112037297B
- Authority
- CN
- China
- Prior art keywords
- grid
- nodes
- system matrix
- tissue
- animal model
- Prior art date
- Legal status (The legal status is an assumption and is not a legal conclusion. Google has not performed a legal analysis and makes no representation as to the accuracy of the status listed.)
- Active
Links
- 230000005540 biological transmission Effects 0.000 title claims abstract description 41
- 238000000034 method Methods 0.000 title claims abstract description 41
- 230000006798 recombination Effects 0.000 title claims abstract description 19
- 238000005215 recombination Methods 0.000 title claims abstract description 19
- 230000003287 optical effect Effects 0.000 claims abstract description 35
- 238000010171 animal model Methods 0.000 claims abstract description 32
- 238000010168 coupling process Methods 0.000 claims abstract description 15
- 230000008878 coupling Effects 0.000 claims abstract description 13
- 238000005859 coupling reaction Methods 0.000 claims abstract description 13
- 239000011159 matrix material Substances 0.000 claims description 87
- 210000001519 tissue Anatomy 0.000 claims description 87
- 230000004907 flux Effects 0.000 claims description 33
- 230000008569 process Effects 0.000 claims description 11
- 238000010521 absorption reaction Methods 0.000 claims description 9
- 210000004185 liver Anatomy 0.000 claims description 6
- 238000009792 diffusion process Methods 0.000 claims description 5
- 230000008521 reorganization Effects 0.000 claims description 5
- 230000015572 biosynthetic process Effects 0.000 claims description 4
- 210000002216 heart Anatomy 0.000 claims description 4
- 210000003734 kidney Anatomy 0.000 claims description 4
- 210000004072 lung Anatomy 0.000 claims description 4
- 210000003205 muscle Anatomy 0.000 claims description 4
- 210000002784 stomach Anatomy 0.000 claims description 4
- 238000003786 synthesis reaction Methods 0.000 claims description 4
- 238000013170 computed tomography imaging Methods 0.000 claims description 3
- 238000012360 testing method Methods 0.000 claims description 3
- 238000003384 imaging method Methods 0.000 abstract description 15
- 238000003325 tomography Methods 0.000 abstract description 12
- 238000004088 simulation Methods 0.000 abstract description 8
- 230000005284 excitation Effects 0.000 abstract description 6
- 238000005415 bioluminescence Methods 0.000 abstract description 4
- 230000029918 bioluminescence Effects 0.000 abstract description 4
- 238000012634 optical imaging Methods 0.000 abstract description 3
- 238000004364 calculation method Methods 0.000 description 7
- 238000010586 diagram Methods 0.000 description 4
- 230000007547 defect Effects 0.000 description 3
- 238000002474 experimental method Methods 0.000 description 3
- 238000013041 optical simulation Methods 0.000 description 3
- 239000000523 sample Substances 0.000 description 3
- 230000005855 radiation Effects 0.000 description 2
- 238000011160 research Methods 0.000 description 2
- 238000005070 sampling Methods 0.000 description 2
- 238000004875 x-ray luminescence Methods 0.000 description 2
- 230000009286 beneficial effect Effects 0.000 description 1
- 230000008859 change Effects 0.000 description 1
- 238000001514 detection method Methods 0.000 description 1
- 238000011161 development Methods 0.000 description 1
- 201000010099 disease Diseases 0.000 description 1
- 208000037265 diseases, disorders, signs and symptoms Diseases 0.000 description 1
- 239000003814 drug Substances 0.000 description 1
- 230000000694 effects Effects 0.000 description 1
- 230000014509 gene expression Effects 0.000 description 1
- 230000003278 mimic effect Effects 0.000 description 1
- 239000003068 molecular probe Substances 0.000 description 1
- 238000010172 mouse model Methods 0.000 description 1
- 230000008520 organization Effects 0.000 description 1
- 230000003252 repetitive effect Effects 0.000 description 1
- 230000035945 sensitivity Effects 0.000 description 1
- 230000005477 standard model Effects 0.000 description 1
- 230000002194 synthesizing effect Effects 0.000 description 1
- 230000001225 therapeutic effect Effects 0.000 description 1
Classifications
-
- G—PHYSICS
- G06—COMPUTING; CALCULATING OR COUNTING
- G06T—IMAGE DATA PROCESSING OR GENERATION, IN GENERAL
- G06T11/00—2D [Two Dimensional] image generation
- G06T11/003—Reconstruction from projections, e.g. tomography
-
- G—PHYSICS
- G06—COMPUTING; CALCULATING OR COUNTING
- G06T—IMAGE DATA PROCESSING OR GENERATION, IN GENERAL
- G06T11/00—2D [Two Dimensional] image generation
- G06T11/003—Reconstruction from projections, e.g. tomography
- G06T11/006—Inverse problem, transformation from projection-space into object-space, e.g. transform methods, back-projection, algebraic methods
-
- G—PHYSICS
- G06—COMPUTING; CALCULATING OR COUNTING
- G06T—IMAGE DATA PROCESSING OR GENERATION, IN GENERAL
- G06T17/00—Three dimensional [3D] modelling, e.g. data description of 3D objects
- G06T17/10—Constructive solid geometry [CSG] using solid primitives, e.g. cylinders, cubes
-
- G—PHYSICS
- G06—COMPUTING; CALCULATING OR COUNTING
- G06T—IMAGE DATA PROCESSING OR GENERATION, IN GENERAL
- G06T17/00—Three dimensional [3D] modelling, e.g. data description of 3D objects
- G06T17/20—Finite element generation, e.g. wire-frame surface description, tesselation
-
- G—PHYSICS
- G06—COMPUTING; CALCULATING OR COUNTING
- G06T—IMAGE DATA PROCESSING OR GENERATION, IN GENERAL
- G06T2210/00—Indexing scheme for image generation or computer graphics
- G06T2210/41—Medical
Landscapes
- Physics & Mathematics (AREA)
- Engineering & Computer Science (AREA)
- General Physics & Mathematics (AREA)
- Theoretical Computer Science (AREA)
- Geometry (AREA)
- Computer Graphics (AREA)
- Software Systems (AREA)
- Mathematical Analysis (AREA)
- Mathematical Optimization (AREA)
- Mathematical Physics (AREA)
- Pure & Applied Mathematics (AREA)
- Algebra (AREA)
- Investigating Or Analysing Materials By Optical Means (AREA)
Abstract
本发明公开了一种基于网格重组的混合光传输建模方法,根据生物组织的光学特性,将各组织划分为高散射组织区域和低散射组织区域;根据各组织所属区域将动物模型离散得到的由空间位置排序的节点和四面体重新排序,得到两个独立网格;基于光传输理论和界面耦合,把两个网格合成一个重组网格,建立混合光传输模型。光传输模型用于前向仿真,能够由光源信息和生物组织结构信息获得生物表面光分布,用于逆向仿真,能够由生物表面光分布和组织结构信息重建光源。本发明能够在精度和效率之间取得更好的平衡,适用于常见光学分子成像模态,包括生物发光断层成像、激发荧光断层成像,契伦科夫光学成像等,在医学分子影像、重建方法等领域有重要的应用价值。
Description
技术领域
本发明属于分子影像领域,涉及一种基于网格重组的混合光传输建模方法,适用于常见的光学分子成像模态中的光传输建模,包括生物发光断层成像、激发荧光断层成像,契伦科夫光学成像、X射线发光断层成像等。
背景技术
光学分子影像是一个快速发展的生物医学成像领域,将光学过程与一定的分子性质相结合,利用特异性探针的光学性质,反映活体状态下分子水平变化,用于组织病理变化的早期检测。该成像方式具有无辐射、反馈快、灵敏度高、时空分辨率高、特异性高等优点,常被用于小动物研究中的分子、细胞和基因表达成像,促进药物开发、疾病研究和治疗干预。
常见的光学分子成像模态,根据成像过程是否需要外部能量激发,分为直接光学分子成像和间接光学分子成像。直接光学分子成像是无需外部能量激发的成像技术,主要通过探测特异性分子探针或放射性核素释放的光子实现成像,包括生物发光断层成像和契伦科夫光学成像。间接光学分子成像是需要外部能量激发的成像技术,如用激光作为外部能量源的荧光分子断层成像,用X射线作为激发源的X射线发光断层成像。
对于各模态光学分子断层成像中,光传输模型用于描述光在生物组织中的传输过程,其准确性和运算效率对光源的快速精确重建起着决定性作用。辐射传输方程(RTE)虽已经被较为广泛地用作描述可见光和近红外光子在生物组织中传播的标准模型。然而在实际应用中,对于复杂的生物组织,RTE的实现是极其复杂的,且需要大量的计算时间。一些RTE近似模型先后被提出,用于模拟光在非均匀介质中的传输,如扩散方程(DE)、简化球谐近似方程(SPN)、离散序数方程(SN)、球谐方程(PN)和相位近似(PA)。然而,DE计算量小、求解速度快,但只适用于高散射特性的生物组织。高阶近似SN、PN和PA,可以精确地模拟光在生物组织内的传输,但较低的计算效率限制了它们在复杂非均匀性介质中的应用。SPN在所有RTE高阶近似模型中,计算量相对较小,更适用于描述光在低散射特性介质中的传输过程。已有的基于扩散方程-简化球谐近似构建的混合光传输模型,求解过程较为复杂,计算量大。光传输模型计算精度和效率的缺陷,制约了光学分子断层成像中光源的快速精确重建。因此,目前需要一种新的光传输模型来实现计算精度和效率之间的平衡。
发明内容
针对现有技术中的缺陷和不足,本发明提供了一种基于网格重组的混合光传输建模方法,解决现有光传输模型模拟光在非均匀介质中的传输时,计算精度和效率之间难以达到有效平衡的问题。
为达到上述目的,本发明采取如下的技术方案:
一种基于网格重组的混合光传输建模方法,该方法根据动物组织的光学特性,将动物各组织划分为高散射组织区域和低散射组织区域;根据各组织所属区域将动物模型离散得到由空间位置排序的节点和四面体重新排序,得到两个独立网格;在这两个网格上,基于光传输理论,高散射组织区域采用扩散方程DE近似建模,低散射组织区域采用三阶简化球谐SP3近似建模,两个区域有相应的相关系统矩阵;通过耦合这两个系统矩阵,把两个网格合成一个重组网格,建立混合光传输模型。
本发明还包括如下技术特征:
具体的,包括以下步骤:
步骤一,获取动物模型并进行预处理,得到动物模型各组织和初始离散网格;
步骤二,根据动物模型各组织的光学参数,并将各组织分为高散射组织区域和低散射组织区域,其中高散射组织区域为Region1,低散射组织区域为Region2;
步骤三,将步骤一得到的初始离散网格中的节点和四面体重新排序,根据节点和四面体所属区域,属于Region1的节点和四面体组成网格1,其中包含的节点数表示为N1;属于Region2的节点和四面体组成网格2,其中包含的节点数表示为N2;
步骤四,在步骤三得到的网格1上,其所属区域作为高散射组织区域采用DE近似建模得到对应的系统矩阵M1;在网格2上,其所属区域作为低散射组织区域采用SP3近似建模得到对应的系统矩阵M2;
步骤五,步骤四得到的系统矩阵M1,维数为N1×N1,N1为网格1包含的节点数;系统矩阵M2,维数为2N2×2N2,N2为网格2包含的节点数;
步骤六,将系统矩阵M1和系统矩阵M2联结组成混合系统矩阵M’,并对重复项耦合后得到耦合后的混合系统矩阵M;
步骤七,把网格1和网格2合成一个整体的重组网格,建立重组网格表面光子通量密度和光源功率密度的关系:
其中,M为步骤六获得的耦合后的混合系统矩阵,Φ为表面光子通量密度,S为光源功率密度,F是源权重矩阵;Φ1、F1和S1对应于网格1;Φ2、F2和S2对应于网格2;合成一个整体的重组网格后,上式表示为:
MΦ=FS
经变换为
Φ=M-1FS
得到混合光传输模型;经过运算可由已知光源信息,获得动物模型表面光子通量密度分布;或已知动物模型表面光子通量密度分布信息,逆向重建出光源信息。
具体的,所述步骤一包括:
步骤1.1,由试验动物的CT成像切片数据中提取组织信息,获得动物模型,该动物模型包括肌肉、肺、心脏、肝脏、胃和肾脏各组织和光源;
步骤1.2,采用Amira软件对动物模型进行剖分,将动物模型离散获得初始离散网格,其中包括N个网格节点、X个内部四面体,该初始离散网格上,节点和四面体是按照其空间位置进行排列。
具体的,所述步骤二包括:
步骤2.1,根据光源发射波长获得各组织对应的光学特性参数,所述光学特性参数包括吸收系数μa,散射系数μs,各向异性系数g和相对折射率n,通过μs′=μs(1-g)计算得到约化散射系数μs′;
步骤2.2,根据光学特性参数将各组织分组为两个区域:高散射组织区域Region1和低散射组织区域Region2,其中,当约化散射系数与吸收系数的比值大于给定阈值时,将该组织作为高散射组织区域;反之,该组织将作为低散射组织区域。
具体的,所述步骤六中,将系统矩阵M1和系统矩阵M2组成一个混合系统矩阵M’:
其中,由于网格1和网格2在边界上相连,边界上的光子通量密度是连续的,即
Φ1′=Φ2′
Φ1’为Region1边界上的光子通量密度,Φ2’为Region2边界上的光子通量密度;
边界上的节点同时属于Region1和Region2,对应于空间同一点,在构成混合系统矩阵时产生了重复项,为避免边界节点的重复贡献,需对其进行耦合,包括以下步骤:
步骤6.1,首先找出网格1和网格2中,空间位置相同的节点,即边界上的重复节点;
步骤6.2,由上述重复节点,获得混合系统矩阵M’中来自于M1和M2部分对应的重复项,即M1部分的一行对应M2部分的两行;
步骤6.3,对重复项进行操作,将M2中重复内容移动到M1对应的行,然后将M2对应项设为0;将整行的项设置为0时,这些共有节点与其自身的关系也被删除,为了恢复这种关系,将上述被设置为0的项中属于系统矩阵M’中对角线上的项设为1;在混合系统矩阵M’中,每个共有节点与其他共有节点的关系表示为:
其中,a对应Region1部分重复节点的项,b和c对应Region2部分重复节点的项;
通过上式对混合系统矩阵M’中对应项进行操作,完成对两个区域连接边界上重复节点的耦合,此过程中,M’中对应项发生变化,可由分块矩阵M11、M12、M21、M22表示,重复项耦合后得到耦合后的混合系统矩阵M,整个过程表示为:
其中,M1为Region1对应的系统矩阵,M2为Region2对应的系统矩阵,M’为混合系统矩阵,M为耦合后的混合系统矩阵;M11、M12、M21、M22为分块矩阵。
本发明与现有技术相比,有益的技术效果是:
本发明通过网格重组和系统矩阵的联结,使混合光传输模型能够一步求解,联结后混合系统矩阵M的维数为N1+2N2;与传统混合模型中,不同的近似方程需要分步求解,混合系统矩阵的维数为有限元体网格总节点数的两倍相比,本发明求解步骤简单,得到的混合系统矩阵维数较低,运算量相对较小,计算效率得到提高。
根据高散射组织区域和低散射组织区域分别采用DE近似模型和SP3近似模型;为每个组织选择合适的光传输模型,提高了计算精度。
附图说明
图1为本发明总体流程图;
图2为网格重组过程的示意图;
图3为一个高度简化的包含2对共有节点的混合网格例子;
图4各模型仿真的表面光子通量密度通过Tecplot软件的展示图;其中图(a)为分子光学模拟环境(MOSE)平台仿真的结果;图(b)为DE近似模型的运行结果;图(c)为SP3近似模型的运行结果;图(d)为本发明提出的混合模型的运行结果。
图5数字鼠实验中各个光传输模型获得的表面光子通量密度比较结果;其中图(a)为各采样点对应光子通量密度的比较;图(b)为各模型对应的平均光子通量密度误差的比较。
具体实施方式
以下对本发明的具体实施方式进行详细说明。应当理解的是,此处所描述的具体实施方式仅用于说明和解释本发明,并不用于限制本发明。
本发明提出了一种基于网格重组的混合光传输建模方法,该方法根据动物组织的光学特性,将动物各组织划分为高散射组织区域和低散射组织区域;根据各组织所属区域将动物模型离散得到由空间位置排序的节点和四面体重新排序,得到两个独立网格;在这两个网格上,基于光传输理论,高散射组织区域采用扩散方程DE近似建模,低散射组织区域采用三阶简化球谐SP3近似建模,两个区域有相应的相关系统矩阵;通过耦合这两个系统矩阵,把两个网格合成一个重组网格,建立混合光传输模型。
包括以下步骤:
步骤一,获取动物模型并进行预处理,得到动物模型各组织和初始离散网格;
步骤1.1,由试验动物的CT成像切片数据中提取组织信息,获得动物模型,该动物模型包括肌肉、肺、心脏、肝脏、胃和肾脏各组织和光源;
步骤1.2,采用Amira软件对动物模型进行剖分,将动物模型离散获得初始离散网格,其中包括N个网格节点、X个内部四面体,该初始离散网格上,节点和四面体是按照其空间位置进行排列。
步骤二,根据动物模型各组织的光学参数,并将各组织分为高散射组织区域和低散射组织区域,其中高散射组织区域为Region1,低散射组织区域为Region2;
步骤2.1,根据光源发射波长获得各组织对应的光学特性参数,所述光学特性参数包括吸收系数μa,散射系数μs,各向异性系数g和相对折射率n,通过μs′=μs(1-g)计算得到约化散射系数μs′;
步骤2.2,根据光学特性参数将各组织分组为两个区域:高散射组织区域Region1和低散射组织区域Region2,其中,当约化散射系数与吸收系数的比值大于给定阈值时,将该组织作为高散射组织区域;反之,该组织将作为低散射组织区域。
步骤三,将步骤一得到的初始离散网格中的节点和四面体重新排序,根据节点和四面体所属区域,属于Region1的节点和四面体组成网格1,其中包含的节点数表示为N1;属于Region2的节点和四面体组成网格2,其中包含的节点数表示为N2;
步骤四,在步骤三得到的网格1上,其所属区域作为高散射组织区域采用DE近似建模;在网格2上,其所属区域作为低散射组织区域采用SP3近似建模;
步骤四中,DE近似建模为:
其中,为梯度,Φ(r)是光密度,S(r)是光源的功率密度,μa为吸收系数,μs为散射系数,g为各向异性系数,An为边界/>处的折射失配系数;ν(r)是边界上的外向单位法向量;
所述SP3近似建模为:
其中为光子通量合成矩,μai=μa+μs(1-gi)(i=1,2,3)为第i阶吸收系数,S(r)是光源的功率密度,μa为吸收系数,μs为散射系数,g为各向异性系数,ν(r)是边界上的外向单位法向量,Ai,Bi,Ci,Di(i=1,2),Ji(i=0,1,2,3)都为系数;
步骤五,由DE方程,获得Region1对应的系统矩阵M1,其维数为N1×N1;由SP3方程,获得Region2对应的系统矩阵M2,其维数为2N2×2N2;
步骤六,将系统矩阵M1和系统矩阵M2联结组成混合系统矩阵M’,并对重复项耦合后得到耦合后的混合系统矩阵M;
具体的,步骤六中,将系统矩阵M1和系统矩阵M2组成混合系统矩阵M’:
其中,由于网格1和网格2在边界上相连,边界上的光子通量密度是连续的,即
Φ1′=Φ2′
Φ1’为Region1边界上的光子通量密度,Φ2’为Region2边界上的光子通量密度;
边界上的节点同时属于Region1和Region2,对应于空间同一点,在构成混合系统矩阵时产生了重复项,为避免边界节点的重复贡献,需对其进行耦合,包括以下步骤:
步骤6.1,首先找出网格1和网格2中,空间位置相同的节点,即边界上的重复节点;
步骤6.2,由上述重复节点,获得系统矩阵M’中来自于M1和M2部分对应的重复项,即M1部分的一行对应M2部分的两行;
步骤6.3,对重复项进行操作,将M2中重复内容移动到M1对应的行,然后将M2对应项设为0;将整行的项设置为0时,这些共有节点与其自身的关系也被删除,为了恢复这种关系,将上述被设置为0的项中属于混合系统矩阵M’中对角线上的项设为1;在混合系统矩阵M’中,每个共有节点与其他共有节点的关系表示为:
其中,a对应Region1部分重复节点的项,b和c对应Region2部分重复节点的项;
通过上式对混合系统矩阵M’中对应项进行操作,完成对两个区域连接边界上重复节点的耦合,此过程中,M’中对应项发生变化,可由分块矩阵M11、M12、M21、M22表示,重复项耦合后得到耦合后的混合系统矩阵M,整个过程表示为:
其中,M1为Region1对应的系统矩阵,M2为Region2对应的系统矩阵,M’为混合系统矩阵,M为耦合后的混合系统矩阵;M11、M12、M21、M22为分块矩阵。
步骤七,把网格1和网格2合成一个整体的重组网格,建立重组网格表面光子通量密度和光源功率密度的关系:
其中,M为步骤六获得的混合系统矩阵,Φ为表面光子通量密度,S为光源功率密度,F是源权重矩阵;Φ1、F1和S1对应于网格1;Φ2、F2和S2对应于网格2;合成一个整体的重组网格后,上式表示为:
MΦ=FS
经变换为
Φ=M-1FS
得到混合光传输模型;经过运算可由已知光源信息,获得动物模型表面光子通量密度分布;或已知动物模型表面光子通量密度分布信息,逆向重建出光源信息。
实施例1:
本实施例采用常用的数字小鼠模型,仅选取小鼠躯干部分,高度35mm,作为研究区域。选用波长为610nm作为本发明中生物发光断层成像的发射波长。光源设置为一个半径为0.8mm、高为1.6mm的圆柱形,放置在肝脏中,中心点坐标(19mm;7mm;17.5mm),光源强度设为1000nw。
将小鼠躯干模型肌肉、心脏、肝脏、肺、胃、肾脏、光源等7个组织,离散为12967个网格节点、67102个内部四面体。并根据光学特性参数分组为高散射特性区域和低散射特性区域,其中肝脏和光源作为低散射特性区域(Region2),其他组织作为高散射特性区域(Region1),如附图1所示。
Region1对应的网格1包含11665个网格节点、55275个内部四面体,作为高散射特性区域采用扩散方程近似(DE)建模,系统矩阵M1的维数为11665×11665。Region2对应的网格2包含2776个网格节点、11827个内部四面体。作为低散射特性区域采用三阶简化球谐近似(SP3)建模,系统矩阵M2的维数为5552×5552。通过对节点空间位置的比较得到两个网格中包含1474个重复节点,对重复节点耦合,获得混合系统矩阵M,维数为17217×17217。把两个网格合成一个重组网格,获得混合光传输模型。
图2为网格重组过程的示意图;其中图(a)为动物模型离散得到的初始网格,包含18个节点,按照其空间位置进行排列;图(b)和图(c)分别为重组后的网格1和网格2,节点根据所属区域重新排序,网格1包含14个节点,网格2包含10个节点。在网格1和网格2的边界上,同属于两个网格的为6个重复节点。
图3为一个高度简化的包含2对共有节点的混合网格例子;进一步解释混合系统矩阵中重复节点对应项的耦合过程:以一组基本网格为例,节点3、4与节点5、6为重复节点;所得混合系统矩阵,用对应的节点编号标记每一行,如附图1中的界面耦合过程;例如,节点3和节点5对应于空间中的同一点,应当被耦合;将混合系统矩阵中节点5的项移动到第3行,然后将节点5中的每一项设为零;为了恢复共有节点与其自身的关系,将节点5中被设为零的项属于混合系统矩阵对角线元素的项设置为1;为了恢复每个共有节点与其他共有节点的关系,节点5的项中连接b和a的元素设置为1,连接b和c的元素设置为2/3。
以分子光学模拟环境(MOSE)平台的仿真结果作为标准,采用单一的DE近似模型和SP3近似模型与本发明提出的基于网格重组的混合光传输模型(HSDM)进行对比。分别采用不同近似建立的模型进行前向仿真,附图4为采用MOSE、DE近似模型、SP3近似模型和HSDM计算的表面光子通量通过Tecplot软件的展示图,其中图(a)为分子光学模拟环境(MOSE)平台仿真的结果;图(b)为DE近似模型的运行结果;图(c)为SP3近似模型的运行结果;图(d)为本发明提出的混合模型的运行结果。由于数值方法计算的结果和MOSE平台仿真的结果单位量纲不同,为了更明确地说明比较,所有结果都进行归一化处理。
为了进一步评估实验结果,附图5展示了实验中各个光传输模型获得的表面光子通量密度在x=19mm处各采样点的结果比较,并计算DE近似模型、SP3近似模型和HSDM的计算结果与MOSE仿真结果的平均相对误差ARE:
其中Simulationi为各模型通过MATLAB仿真得到的表面光子通量密度,Mosei为通过MOSE平台得到的表面光子通量密度,N为采样点总数。
图5(a)为各采样点对应光子通量密度的比较;图(b)为各模型对应的平均光子通量密度误差的比较,从附图5可以明显看出,与DE近似模型和SP3近似模型相比,HSDM的结果更接近于MOSE的结果。DE、SP3和HSDM模型的ARE分别为0.3011,0.2479和0.1606,HSDM的平均相对误差最小。
本实施例中,数字小鼠实验得到了表面光子通量密度,同时记录了各模型运行消耗的时间。DE近似模型平均用时103.57s,SP3近似模型平均用时561.89s,HSDM平均用时247.81s,计算比较得HSDM比SP3近似模型用时减少55.9%。
实验结果表明,与DE近似模型和SP3近似模型相比,本发明提出的基于网格重组的混合光传输模型在精度和效率之间取得了更好的平衡,相对更适合作为光学分子断层成像的光传输模型。
Claims (3)
1.一种基于网格重组的混合光传输建模方法,其特征在于,该方法根据动物组织的光学特性,将动物各组织划分为高散射组织区域和低散射组织区域;根据各组织所属区域将动物模型离散得到由空间位置排序的节点和四面体重新排序,得到两个独立网格;在这两个网格上,基于光传输理论,高散射组织区域采用扩散方程DE近似建模,低散射组织区域采用三阶简化球谐SP3近似建模,两个区域有相应的相关系统矩阵;通过耦合这两个系统矩阵,把两个网格合成一个重组网格,建立混合光传输模型;
包括以下步骤:
步骤一,获取动物模型并进行预处理,得到动物模型各组织和初始离散网格;
步骤二,根据动物模型各组织的光学参数,并将各组织分为高散射组织区域和低散射组织区域,其中高散射组织区域为Region1,低散射组织区域为Region2;
步骤三,将步骤一得到的初始离散网格中的节点和四面体重新排序,根据节点和四面体所属区域,属于Region1的节点和四面体组成网格1,其中包含的节点数表示为N1;属于Region2的节点和四面体组成网格2,其中包含的节点数表示为N2;
步骤四,在步骤三得到的网格1上,其所属区域作为高散射组织区域采用DE近似建模得到对应的系统矩阵M1;在网格2上,其所属区域作为低散射组织区域采用SP3近似建模得到对应的系统矩阵M2;
步骤五,步骤四得到的系统矩阵M1,维数为N1×N1,N1为网格1包含的节点数;系统矩阵M2,维数为2N2×2N2,N2为网格2包含的节点数;
步骤六,将系统矩阵M1和系统矩阵M2联结组成混合系统矩阵M’,并对重复项耦合后得到耦合后的混合系统矩阵M;
步骤七,把网格1和网格2合成一个整体的重组网格,建立重组网格表面光子通量密度和光源功率密度的关系:
其中,M为步骤六获得的耦合后的混合系统矩阵,Φ为表面光子通量密度,S为光源功率密度,F是源权重矩阵;Φ1、F1和S1对应于网格1;
Φ2、F2和S2对应于网格2;合成一个整体的重组网格后,上式表示为:
MΦ=FS
经变换为
Φ=M-1FS
得到混合光传输模型;经过运算可由已知光源信息,获得动物模型表面光子通量密度分布;或已知动物模型表面光子通量密度分布信息,逆向重建出光源信息;
所述步骤六中,将系统矩阵M1和系统矩阵M2组成一个混合系统矩阵M’:
其中,由于网格1和网格2在边界上相连,边界上的光子通量密度是连续的,即
Φ1′=Φ2′
Φ1’为Region1边界上的光子通量密度,Φ2’为Region2边界上的光子通量密度;
边界上的节点同时属于Region1和Region2,对应于空间同一点,在构成混合系统矩阵时产生了重复项,为避免边界节点的重复贡献,需对其进行耦合,包括以下步骤:
步骤6.1,首先找出网格1和网格2中,空间位置相同的节点,即边界上的重复节点;
步骤6.2,由上述重复节点,获得混合系统矩阵M’中来自于M1和M2部分对应的重复项,即M1部分的一行对应M2部分的两行;
步骤6.3,对重复项进行操作,将M2中重复内容移动到M1对应的行,然后将M2对应项设为0;将整行的项设置为0时,这些共有节点与其自身的关系也被删除,为了恢复这种关系,将上述被设置为0的项中属于系统矩阵M’中对角线上的项设为1;在混合系统矩阵M’中,每个共有节点与其他共有节点的关系表示为:
其中,a对应Region1部分重复节点的项,b和c对应Region2部分重复节点的项;
通过上式对混合系统矩阵M’中对应项进行操作,完成对两个区域连接边界上重复节点的耦合,此过程中,M’中对应项发生变化,可由分块矩阵M11、M12、M21、M22表示,重复项耦合后得到耦合后的混合系统矩阵M,整个过程表示为:
其中,M1为Region1对应的系统矩阵,M2为Region2对应的系统矩阵,M’为混合系统矩阵,M为耦合后的混合系统矩阵;M11、M12、M21、M22为分块矩阵。
2.如权利要求1所述的基于网格重组的混合光传输建模方法,其特征在于,所述步骤一包括:
步骤1.1,由试验动物的CT成像切片数据中提取组织信息,获得动物模型,该动物模型包括肌肉、肺、心脏、肝脏、胃和肾脏各组织和光源;
步骤1.2,采用Amira软件对动物模型进行剖分,将动物模型离散获得初始离散网格,其中包括N个网格节点、X个内部四面体,该初始离散网格上,节点和四面体是按照其空间位置进行排列。
3.如权利要求1所述的基于网格重组的混合光传输建模方法,其特征在于,所述步骤二包括:
步骤2.1,根据光源发射波长获得各组织对应的光学特性参数,所述光学特性参数包括吸收系数μa,散射系数μs,各向异性系数g和相对折射率n,通过μs′=μs(1-g)计算得到约化散射系数μs′;
步骤2.2,根据光学特性参数将各组织分组为两个区域:高散射组织区域Region1和低散射组织区域Region2,其中,当约化散射系数与吸收系数的比值大于给定阈值时,将该组织作为高散射组织区域;反之,该组织将作为低散射组织区域。
Priority Applications (1)
Application Number | Priority Date | Filing Date | Title |
---|---|---|---|
CN202010742538.5A CN112037297B (zh) | 2020-07-29 | 2020-07-29 | 一种基于网格重组的混合光传输建模方法 |
Applications Claiming Priority (1)
Application Number | Priority Date | Filing Date | Title |
---|---|---|---|
CN202010742538.5A CN112037297B (zh) | 2020-07-29 | 2020-07-29 | 一种基于网格重组的混合光传输建模方法 |
Publications (2)
Publication Number | Publication Date |
---|---|
CN112037297A CN112037297A (zh) | 2020-12-04 |
CN112037297B true CN112037297B (zh) | 2023-09-22 |
Family
ID=73583452
Family Applications (1)
Application Number | Title | Priority Date | Filing Date |
---|---|---|---|
CN202010742538.5A Active CN112037297B (zh) | 2020-07-29 | 2020-07-29 | 一种基于网格重组的混合光传输建模方法 |
Country Status (1)
Country | Link |
---|---|
CN (1) | CN112037297B (zh) |
Citations (2)
Publication number | Priority date | Publication date | Assignee | Title |
---|---|---|---|---|
WO2009044282A2 (en) * | 2007-10-04 | 2009-04-09 | Mental Images Gmbh | Quasi-monte carlo light transport simulation by efficient ray tracing |
CN102393969A (zh) * | 2011-06-02 | 2012-03-28 | 西安电子科技大学 | 基于生物组织特异性的光学三维成像方法 |
Family Cites Families (1)
Publication number | Priority date | Publication date | Assignee | Title |
---|---|---|---|---|
CN104921706B (zh) * | 2015-07-08 | 2018-02-09 | 北京工业大学 | 基于多任务贝叶斯压缩感知方法的生物发光断层成像重建算法 |
-
2020
- 2020-07-29 CN CN202010742538.5A patent/CN112037297B/zh active Active
Patent Citations (2)
Publication number | Priority date | Publication date | Assignee | Title |
---|---|---|---|---|
WO2009044282A2 (en) * | 2007-10-04 | 2009-04-09 | Mental Images Gmbh | Quasi-monte carlo light transport simulation by efficient ray tracing |
CN102393969A (zh) * | 2011-06-02 | 2012-03-28 | 西安电子科技大学 | 基于生物组织特异性的光学三维成像方法 |
Non-Patent Citations (1)
Title |
---|
任晓芳 ; 徐亮 ; .利用辐射传输-扩散混合模型的DOT图像重建算法.微型电脑应用.2016,(02),全文. * |
Also Published As
Publication number | Publication date |
---|---|
CN112037297A (zh) | 2020-12-04 |
Similar Documents
Publication | Publication Date | Title |
---|---|---|
CN103271723B (zh) | 一种生物发光断层成像重建方法 | |
Joshi et al. | Radiative transport-based frequency-domain fluorescence tomography | |
CN108451508B (zh) | 基于多层感知机的生物自发荧光三维成像方法 | |
CN111915733B (zh) | 基于LeNet网络的三维锥束X射线发光断层成像方法 | |
Cong et al. | Differential evolution approach for regularized bioluminescence tomography | |
Cong et al. | A Born‐type approximation method for bioluminescence tomography | |
Wang et al. | A novel adaptive parameter search elastic net method for fluorescent molecular tomography | |
WO2017004851A1 (zh) | 基于多任务贝叶斯压缩感知方法的生物发光断层成像重建算法 | |
CN110680284A (zh) | 一种基于3D-Unet的介观荧光分子成像三维重建方法及系统 | |
CN1916923A (zh) | 人体内物质代谢功能信息可视化的计算机模拟方法 | |
CN109166103A (zh) | 基于多层感知网络的激发荧光断层成像方法 | |
Cong et al. | A practical method to determine the light source distribution in bioluminescent imaging | |
CN112037297B (zh) | 一种基于网格重组的混合光传输建模方法 | |
Feng et al. | An adaptive regularization parameter choice strategy for multispectral bioluminescence tomography | |
Jiang et al. | Image reconstruction for bioluminescence tomography | |
CN115778412A (zh) | X光声成像中造影剂剂量的优化方法、装置及存储介质 | |
Mathews et al. | A generalized reconstruction framework for unconventional PET systems | |
Kumar et al. | Monte Carlo method for bioluminescence tomography | |
CN107374588B (zh) | 一种基于同步聚类的多光源荧光分子断层成像重建方法 | |
El Fakhri et al. | A new scatter compensation method for Ga-67 imaging using artificial neural networks | |
Zhang et al. | A meshfree method for solving cardiac electrical propagation | |
CN117830373B (zh) | 基于加速预条件邻近梯度算法的pet图像重建方法 | |
CN108846790A (zh) | 一种加速图像重建的方法 | |
Cao et al. | FISTA-NET: Deep Algorithm Unrolling for Cerenkov luminescence tomography | |
Cong et al. | Reconstruction of optical parameters for molecular tomographic imaging |
Legal Events
Date | Code | Title | Description |
---|---|---|---|
PB01 | Publication | ||
PB01 | Publication | ||
SE01 | Entry into force of request for substantive examination | ||
SE01 | Entry into force of request for substantive examination | ||
GR01 | Patent grant | ||
GR01 | Patent grant |