专利名称：声速不均匀介质热声成像的重建算法的制作方法
技术领域：
本发明属于热声成像技术领域，具体为一种声速不均匀介质热声成像的图像重建算法。技术背景热声成像(往往也称为光声成像)结合了光学成像高对比度和超声成像高分辨率的优 点，近年来发展迅速，成为一种广泛关注的新型无损医学成像技术。热声成像非常适用于 肿瘤检测、血管成像等方向，也被尝试应用于老鼠脑部的结构和功能成像、分子成像、流 速检测等方向。热声成像中，生物组织受电磁波(常用微波或激光)脉冲照射后产生热膨胀，从而激发 出热声波。热声波是一种超声波，由超声换能器探测后，可用来重建生物组织内部的电磁 波吸收分布图像。热声成像技术的关键是成像方法，即如何由探测到的热声波重建出组织 的电磁波吸收分布。目前已提出很多成像的重建算法，但是这些算法几乎都假定生物组织 内部的声速是均匀的。事实上，这个假设有时是不成立的，例如对乳房成像时声速差异可达10%。此时若假定声速均匀来进行图像重建，将会导致图像内目标的错位和模糊[1'2]。近两年已提出补偿声速不均匀性的热声成像算法[1—4。这些算法大都需要声速分布的先 验知识[1—3]，而实际应用中往往是不知道声速分布的。目前仅有的无需声速分布先验知识的方法[41是基于有限元法、通过迭代方法进行计算，重建速度很慢。本发明不仅无需介质内声速分布的先验知识，而且是一种直接的反向方法，可以一次性快速重建出声速不均匀介质的整幅热声图像。发明内容本发明的目的在于提出一种无需介质内声速颁布先验知识，速度快的声速不均匀介质热声成像的重建算法。本发明提出的一种声速不均匀介质热声成像的重建方法。其具体步骤为对声速均匀 介质热声成像的时域重建法进行修正，使声波在介质中两点间的传播时间不再与距离成正 比，而是一个待定值；估计出介质中的声速分布，通过微扰理论计算出声波在介质中两点 间的传播时间，然后进行图像重建。下面对各歩骤作进一歩具体描述。通常电磁波脉冲持续时间远短于组织的热扩散时间，因此可忽略热扩散。下面给出三维下热声成像的基本方程[2'5]:<formula>formula see original document page 4</formula>其中/ (r，,)是位置r处的声压，^(r)是电磁波吸收分布，/是时间，/(,)是电磁波能量脉冲函 数，^是等压膨胀系数，"r)是声速分布，Cp是样品组织的比热。热声成像算法是一个典型 的逆问题，即如何由p(r力求出J(r)。当超声换能器沿圆球面轨迹ro(球心为原点，半径为^)扫描时，p(r, 0中r = rQ。通常 假定电磁波脉冲函数是一个^函数，即/(0 = 3(/);声速是一个常数，即c(r)-c。此时，热 声成像的时域重建法为[5]:<formula>formula see original document page 4</formula>其中7是一个常数，O)是积分面。事实上，(2)式的物理意义非常明确。位置r处产生的热声波到达超声换能器位置r0 需要的时间为/。因此X(r)等于所有探测位置的热声波在对应时间f上的值(l")(^(r。力/30进行投影。如果组织中声速是均匀的，f和r与ro间的距离成正比。如果声速不是均匀的， /是关于r与rQ的函数，即/ = T(r， rG)。因此(2)式被修正为<formula>formula see original document page 4</formula>3)(3)式就是本发明的反向重建公式，它考虑到了热声成像中的声速不均匀性。那么现在的问 题是如何求出r(r,r0)。生物软组织内的声速差异一般不大，因此组织的慢度w(r)- 1/c(r)可以写成w(r) = w0 +6>,w1(r) (4) 其中的是常数，e是一个很小的值。根据微扰理论，r(r,ro)可以写成微扰级数的形式"'"<formula>formula see original document page 4</formula> (5) 此处忽略了 二阶以上的小量。可以进 一 步求出[7 ]:<formula>formula see original document page 4</formula> (6)其中/是r和r。间的线段，即声速均匀时的声波传播路径。从(6)式可见，rQ是没有慢度扰动Wl时的声波传播时间；L是一阶时间扰动，表示了^沿/的影响；『2是二阶时间扰动，表示了声波折射的主要影响。显然r2的计算较复杂，因此在生物软组织的热声成像中，考虑到声速差异较小，可以忽略T2的影响，通过To和 7V计算r就足够了。但当介质中的声速差异较大时，例如对生物体的脑部成像时，就需要考虑r2的影响。(6)式需要知道慢度(或声速)分布w(r)(或c(r)),而实际应用中往往不知道w(r)。因此本 发明提出了一种估计w(r)的方法，从而进一步计算出r(r，ro)。具体步骤如下。 超声换能器接收到的声波可写为[2]:(7)7\r,r0)=,其中W是一个常数。定义.-S(r0,f)=」 nn(8)可见S(r。， 0是J(r)的面积分，积分面O)上所有点到r。的声波传播时间为r。若c(r)是常数， 则积分面Dc可简化为球面。如图1所示，当待测组织远比r和ro间距离小时，O)可近似 为一个平面。那么显然有S(r。，,)"(-r。力 (9) 其中附'=r。)。若声速有较小差异，或者待测组织较大时，(9)仍近似成立。则r(-rQ, r0)可通过K式计算:r(-ro'ro)二argma來。(0) (1°) 其中^。 (0是S(r。， 0和S(-r。， -O的相关函数，定义为7 r。 ( ) = ￡ S(r。, * r。 ，-(r - f麵 (11) 可以认为-ro与r。间的声速是近似均匀的，则有<formula>formula see original document page 5</formula>(12)其中&(r)表示估计的慢度分布u (r)。综上，结合(3)、 (6)、 (10)、 (12)式，即可完成声速不均匀的图像重建。首先根据(IO)、 (12)式用探测到的热声波估计出慢度(或声速)分布；然后根据(6)式计算组织内两点间的声 波传播时间；最后根据(3)式重建出待测组织的图像。实际应用中参数e可随意选取，其结 果都一样；参数;y、 ^应根据待测组织的等压膨胀系数、比热、声速以及扫描半径来确定[5]。
图1、 r和r。位置处接收到的热声波所对应的待测组织部分。图2、待测组织的电磁波吸收分布(a)和声速分布(b)。图3、使用不同声速的图像重建结果，(a) c = 1500 m/s， (b) c = 1429 m/s， (c) c = 1364 m/s, (d)已知声速分布。图4、声速分布未知时，估计的声速分布(a)和重建的电磁波吸收分布(b)。
具体实施方式
为验证本发明的方法，在计算机上进行二维情况下的仿真实验。虽然本发明的理论论 证是在三维情况下进行的，但将其应用到二维下显然也是成立的。仿真实验的具体步骤和每歩的结果如下所述1、 建立待测样品的模型。图2显示了待测组织的电磁波吸收分布和声速分布(m/s)。该声速分布类似于女性乳房模型[1]。图2的尺寸为18mmX18mm。待测组织的其 余参数(如等压膨胀系数、比热)对于成像算法的影响体现在常数V、 ^上，其选取可以不用 与真实情况 -致，此处简单设置为l。2、 利用时域有限差分法仿真出扫描圆周上接收到的声波，扫描圆周半径为9mm, 共有160个等间距的探测位置。3、 先假设声速分布是己知的，用仿真的声波重建出待测样品的图像。图3显示了 用均匀声速和已知不均匀声速分布的图像重建结果。4、 声速分布未知时，进行图像重建。先利用本发明的(10)式和(12)式来估计声速 分布，然后再基于(3)式和(6)式进行图像重建。估计的声速分布和重建的电磁波吸收分布图 像如图4所示。由歩骤3的结果可见，图3(a) (c)使用均匀声速，结果不好。其中以图3(b)效果相对最 好，但最小的圆形组织已经模糊。图3(d)使用本发明的(3)式和(6)式进行图像重建，效果很 好。这说明本发明的修正的时域重建法是有效的。重建时，仅对扫描圆周内的区域进行图 像重建，而扫描圆周外在重建时像素直接设为O。由于扫描圆周内区域在重建时存在误差， 像素本应为O的背景区域重建出的结果是很小的负值。因此图3边缘处有一圆形轮廓，这 仅表明扫描轨迹，并不是反映待测物体的轮廓。山歩骤4的结果可见，图4(a)中估计的声速分布与实际情况差别较大，但用在重建电 磁波吸收分布图像时，可以较准确地计算出(6)式，获得满意的图像重建结果。这与本发明 的理论相符。图4(b)中组织轮廓较清晰，最小的圆形组织也基本没有模糊。总的来看，图 4(b)的结果略差于图3(d),好于图3(a) (c)。这是因为图4(b)是基于估计出的声速分布进行 重建，N定不及已知声速分布的图3(d)，但会优于对声速不均匀性不进行任何补偿的图 3(a) (c)。根据实验结果可见，本发明的方法无需声速分布的先验知识，而且可以一次性快速重 建出声速不均匀介质的整幅热声图像，效果很好。本发明的仿真实验中声速差异约为10%, 一般生物软组织的声速差异会更小一些，因此本方法非常适用于生物软组织的热声成像。 参考文献[1] Y. Xu， and L. V. Wang， "Effects of acoustic heterogeneity in breast thermoacoustic tomography,"7Va肌Wfmy肌,Fmw/ect., />叫Cbn/X, vol. 50， no. 9, Sep. 2003.[2] X. Jin, and L. V. Wang， "Thermoacoustic tomography with correction for acoustic speed variations,"Afoi肠/.， vol. 51， pp. 6437-6448， 2006.[3] J. Zhang， and M. A. Anastasio, "Reconstruction of speed-of-sound and electromagnetic absorption distributions in photoacoustic tomography," in 尸n vol. 6080， pp. 608619-1-7,2006.[4] H. Jiang, Z. Yuan, and X. Gu， "Spatially varying optical and acoustic property reconstruction using finite-element-based photoacoustic tomography," /. 6 / /1. vol. 23， no. 4, pp.878-888, Apr. 2006.[5] M. Xu， Y. Xu, and L. V. Wang， "Time-domain reconstruction algorithms and numericalsimulations for thermoacoustic tomography in various geometries," 5z'o附ed￡"g,， vol. 50， no. 9， pp. 1086-1099， Sep. 2003. [6] R. Snieder, and M. Sambridge, "Ray perturbation theory for traveltimes and ray paths in 3-Dheterogeneous media," G^op一. 《/ /《vol. 109， pp. 294-322, 1992. [7] R. Snieder, and D. F. AIdridge， "Perturbation theory for travel times," / JcowW. J州.，vol. 98, no, 3， pp. 1565-1569, Sep. 1995.
权利要求
1、一种声速不均匀介质热声成像的图像重建方法，其特征在于对声速均匀介质热声成像的时域重建法进行修正，使声波在介质中两点间的传播时间不再与距离成正比，而是一个待定值；估计出待测介质中的声速分布，通过微扰理论计算出声波在介质中两点间的传播时间，然后进行图像重建；其中所述的修正的时域重建法为
全文摘要
本发明属于热声成像技术领域，具体为一种声速不均匀介质热声成像的图像重建方法。该方法通过对声速均匀介质热声成像的时域重建法进行修正，使声波在介质中两点间的传播时间不再与距离成正比，而是通过估计出的声速分布计算出来。当前对声速不均匀介质热声成像的重建算法大都需要声速分布的先验知识，并采用迭代方法计算，速度很慢。本发明无需介质内声速分布的先验知识，而且是一种直接的反向方法，可以一次性快速重建出声速不均匀介质的整幅热声图像。
文档编号A61B8/00GK101214156SQ20081003250
公开日2008年7月9日 申请日期2008年1月10日 优先权日2008年1月10日
发明者弛 张, 汪源源 申请人:复旦大学