X射线荧光CT成像中荧光产额、退激时间、散射、偏振等关键物理问题计算与分析*

张芝振 李亮†

1) (清华大学工程物理系,北京 100084)

2) (清华大学,粒子技术与辐射成像教育部重点实验室,北京 100084)

X 射线荧光CT(X-ray fluorescence computed tomography,XFCT)是一种使用X 射线荧光(X-ray fluorescence,XRF)实现功能性成像的新技术,在生物医学成像中表现出较大潜力.但是,X 射线穿过生物体的同时还会产生大量康普顿散射光子,对XRF 信号的采集形成很强的背景噪声;因此,如何有效消除康普顿散射噪声对于提高XFCT 成像质量至关重要.本文研究总结了XFCT 成像过程中涉及的物理过程,包括:荧光的产额、退激发时间、荧光发射角分布、荧光偏振态、康普顿散射角分布与散射光偏振态,并通过研究荧光与散射光物理性质的差异寻找去除康普顿散射噪声的方法.经过物理过程推导和分析计算,发现:1) 高原子序数元素的K 层荧光退激发时间极短,在现有探测器的时间分辨率条件下,无法分辨散射光与荧光;2) K 层发射荧光的角分布各向同性,康普顿散射角分布在与入射光偏振方向附近取得最小值,而且入射光线偏振度越高,散射光的微分截面越小,偏振光源将有利于减少康普顿散射噪声;3) K 层荧光线偏振度为零,而康普顿散射光子在一些散射方向上具有一定线偏振度,因此偏振态的差异可能用于降低康普顿散射噪声.

1 引言

X 射线荧光CT(X-ray fluorescence computed tomography,XFCT)是一种针对高原子序数元素的功能成像技术.XFCT 通过探测X 射线与物体内高原子序数元素发生光电效应后产生的XRF 得到物质内该元素的浓度分布图像,所以兼具透射型CT 结构成像和发射型CT 功能成像的优点.Boisseau 和Grodzins[1]利用美国国家同步辐射光源首次实现XFCT 实验.同步辐射光源具有准单能和高线偏振度的特性,是XFCT 成像的理想光源,但是因为同步辐射占用场地和成本高昂的问题,所以基于同步辐射的XFCT 系统很难应用于实际医学应用中.X 光机相比于同步辐射具有巨大的灵活性和经济性优势,近年来有很多基于X 光机的XFCT 系统设计的研究[2−9].因为X 光机光源采用电子轫致辐射的方式产生X 射线,所以其能谱是分布范围较大的连续谱并且不具备线偏振特性.相比于同步辐射光源,X 光机光源会产生更大的散射背景,严重影响了XFCT 的成像质量.

去除散射背景的方法主要有两种.第一种方法采用略高于和略低于荧光物质的K-edge 的单能入射光分别采取两组数据,两组数据之差可近似为荧光信号[10].第一种方法对光源单能性要求很高,一般采用同步辐射源.第二种方法利用散射光子能量连续分布的特点,通过三次多项式拟合的方式估计散射,进而在数据中减去散射的贡献[2],该方法对探测器的能量分辨率要求高,一般使用单像素X 射线谱仪采集荧光信号,导致扫描时间很长.

逆康普顿散射源相比于同步辐射源不仅具有占地面积小和价格适中的优势,而且具有同步辐射源准单能和线偏振的特性,适合应用于XFCT 成像系统设计中.逆康普顿散射源与同步辐射源降低XFCT 散射背景的原理都是应用了线偏振X 光康普顿散射的角分布特性.Chi 等[11]蒙卡模拟研究了逆康普顿光源线偏振特性对XFCT 成像质量的影响,发现逆康普顿光源能抑制XFCT 成像的康普顿散射背景约1.6 倍.

此外,荧光与散射光的物理性质的差别可能用于降低散射背景,比如是否有可能利用荧光与散射光飞行时间、角分布和偏振态的差别减少散射背景或寻找散射背景最小的探测方位.本文分析了XFCT中荧光与散射光相关物理过程的性质:荧光与散射光的角分布,荧光退激发时间与荧光产额,荧光与散射光的偏振态.这些物理过程决定了XFCT 荧光成像的极限,通过这些物理过程性质能够指导XFCT 的系统设计,减少XFCT 的康普顿散射背景和提高XRF 信号强度.

2 XFCT 成像方法和物理原理

在XFCT 所涉及的能量范围内,X 射线与物质相互作用方式有光电效应、相干散射和康普顿散射.在X 射线能量为几十keV 时,相干散射截面与康普顿散射截面相当,所以在处理XFCT 散射问题时相干散射和康普顿散射一般都要考虑,其中相干散射主要沿着小角方向散射,对光子的传播会有影响,但通过后面的讨论可以看出,在大角度方向康普顿散射是主要的.

XFCT 成像系统几何设计有很多种方式[2−12],其中基于小孔成像、扇束X 射线源和XRF 探测器阵列的成像系统设计能大幅提高扫描速度[5],其XFCT 实验布置示意图如图1 所示.小孔准直器作用是使XRF 探测器具有空间分辨能力.示意图中被扫描物体是一个包含四个小圆柱的大圆柱模型,四个小圆柱里含有不同浓度的某种高原子序数元素的溶液.图1 中描绘了某一条光路以及该路径上可能产生的XRF 光子和散射光子,入射光与溶液中的重金属元素原子发生光电效应,重金属原子发射荧光(特征X 射线),荧光经过小孔后,被XRF探测器阵列探测,同时散射光子也会进入XRF 探测器,形成散射背景.XRF 探测器可通过设置能窗阈值得到入射到每个探测器单元的相应能窗的光子计数.X 射线探测器阵列用于传统的透射型CT成像,其得到的衰减系数分布图可用于XFCT 成像的衰减校正.

图1 基于小孔成像方式的XFCT 实验设置Fig.1.Experimental setup based on pinhole imaging.

为描述入射和散射X 射线的偏振态,引入Stokes矢量[13]:

其中I表示光强度;P1和P2表示在与传播方向垂直的平面内成 4 5◦夹角的两对正交方向下分别测量得到的平面偏振度;P3表示圆偏振度.光束的线偏振度和圆偏振度定义分别为[14]

其中PL表示线偏振度;PC表示圆偏振度;Imax和Imin分别表示经过理想线偏振片后光强的最大与最小值;IR和IL分别表示光束中右旋光和左旋光的强度.定义归一化Stokes 矢量为

Stokes 矢量与正交方向选择有关,选定ξ1−1 对应的线偏振方向为x轴,则ξ1+1 对应的线偏振方向为y轴,传播方向为z轴正向,并规定:迎着光子传播方向(z轴正向)看时,绕z轴顺时针旋转角度为正,绕z轴逆时针旋转角度为负,则在旋转ψ角度后的新正交方向下的Stokes 矢量 [ 1,ζ]T与原正交方向下的Stokes 矢量 [ 1,ξ]T由变换矩阵M联系,即 [ 1,ζ]TM[1,ξ]T.M的非零元为

一束光经过某个相互作用后,其Stokes 矢量会改变,通过一个 4×4 变换矩阵T联系相互作用前后的Stokes 矢量:

3 K 层荧光物理性质

3.1 荧光产额

某个壳层或子壳层的荧光产额定义为该壳层产生一个空穴后,通过发射荧光退激发的概率.K 壳层产生空穴后,不同跃迁方式产生的荧光的能量和强度不同.一般地,壳层越高的电子通过退激发填充K 层空穴,产生的荧光的能量越高,某种荧光的强度正比于该荧光对应跃迁的分支比.荧光产额的半经验公式[15]:

式中,ωK表示K 层荧光产额;Z表示原子序数;Ci表示拟合参数.(6)式中参数Ci可以由实验得到的不同元素荧光产额数据拟合得到.Bambynek 等[15]和Hubbell 等[16]结合理论计算与实验给出拟合参数和元素荧光产额数据表.

K 层荧光产额ωK是原子序数Z的增函数,对Z >60 的元素,ωK>90%,所以高原子序数原子的K 层空穴态通过发射荧光退激的概率最大.因此XFCT 应该选择高原子序数元素作为荧光元素,在高Z元素范围内荧光产额基本饱和,提高Z对增加荧光产额的收益很低,但提高Z可增大荧光能量,能量高的光子更容易穿透物体被探测器采集到,从而降低统计噪声.

K 壳层空穴态不同跃迁方式产生不同能量的荧光,荧光按末态空穴所在层分成两组,即Kα代表来自L 层电子的跃迁产生的荧光,Kβ表示来自M 层及更高层电子跃迁产生的荧光,同一组内的荧光能量相近.Ertuğral 等[17]给出了原子序数16≤Z≤92 范围内59 种元素的Kβ/Kα荧光强度比测量值,Scofield[18]给出了10≤Z≤98 原子序数范围内Kβ/Kα荧光强度比理论计算值,从实验和理论计算结果可以看出,荧光强度比Kβ/Kα随原子序数增大而增大,对于Z>60 的元素,Kβ/Kα荧光强度比在0.3 左右.通过荧光强度比和K 层荧光产额ωK可以得到K 层不同能量荧光的产额.

入射X 射线也会与重金属原子的L,M 和N等更高电子壳层的电子发生光电效应,且发射的荧光能量随层数增高而减小,其中荧光能量最高的L 层荧光的能量一般为十几keV,属于低能区的光子,几乎无法穿透模体,所以L,M 和N 等更高电子壳层的荧光在XFCT 荧光成像中可以忽略.

3.2 荧光退激发时间

如果荧光退激发半衰期很长,且长于探测器的时间分辨率,那么可以采用极短时间的X 光脉冲作为光源,在散射光子到达探测器后,延迟一段时间采集荧光光子,这样就可以在时间上分开散射和荧光从而去除康普顿散射背景;如果探测器的时间分辨率大于荧光退激发半衰期,那么荧光和散射光子到达探测器的平均时间差小于探测器时间分辨率,从而无法分辨荧光和散射光子.为了证明这种方法是否可行,需要计算荧光退激发时间.原子发生光电效应后处于激发态,激发态是不稳定的,假设t=0 时刻原子处于激发态,则在以后的t时刻体系仍然停留在这个态的概率是

式中,τ是激发态平均存在时间,亦称激发态寿命,其与激发态能级宽度Γ关系为

K层荧光分宽度ΓR与K层荧光产额的关系是[15]

对于Z>40 的元素,K层空穴态能级宽度半经验公式为[15]

通过(8)式和(9)式得到K层荧光寿命τK为

可见,K 层荧光寿命随原子序数增加快速下降,对已知的XFCT 应用的最低原子序数的元素钆(Z64),K 层荧光寿命为0.03 ps.当前对伽马光子探测器时间分辨率要求最高的PET 设备的符合时间分辨率在百皮秒量级[19],远大于中高原子序数元素的K 层荧光退激发时间,因此,无法通过时间差别来分辨荧光与散射光子进而减少散射背景.

3.3 荧光发射角分布与线偏振度

原子内壳层荧光发射角分布和偏振态在理论和实验上已经有很多研究[20−26]:若空穴态的总角动量量子数J1/2,则荧光发射角分布各向同性且线偏振度为零[20];若J>1/2,则荧光发射一般各向异性并且部分偏振[22].因为K 层空穴态J1/2[22],所以K 层空穴态发射荧光的角分布各向同性且线偏振度为零.L,M 等不同子层具有不尽相同的总角动量量子数,而且某些子层的空穴态会通过Coster-Kronig 跃迁将空穴转移到其他子层,因此L,M 等壳层的空穴态发射荧光更加复杂多样,剔除Coster-Kronig 效应的影响后,不同子层空穴态发射荧光角分布与线偏振度也满足上述角动量规则[24].

4 XFCT 中的散射

4.1 散射光角分布

实验室参考系下单个自由电子与线偏振光的康普顿散射截面Klein-Nishina 公式为[25]

式中,r0是经典电子半径;k和k0分别是出射和入射光子能量,以电子静质量能为单位;e'和e分别是出射光偏振方向和入射光偏振方向; dΩ是立体角元.假设在出射方向n放置偏振器,该偏振器透光轴沿e'方向,一个偏振态为e的入射光子经过康普顿散射到出射方向n,则该光子通过偏振器的概率正比于n方向上的 dσKN/dΩ[27].散射方向和偏振方向如图2 所示,入射光沿Z轴正向入射,用矢量k0(|k0|k0)表示;散射光沿方向出射,用矢量k(|k|k)表示.e是描述入射光偏振方向的单位矢量,θ∠SOZ是散射角,φ∠AOX是散射的方位角.y和y'表示垂直于散射平面SOX的线偏振方向,x和x'表示平行于平面SOX的线偏振方向,位于平面x′Sy′内的e'表示散射光的偏振方向,e'与y'夹角为β,显然,x',y'和e'垂直kf.从图2描述的几何关系得到 (e·e′)2cos2β(1−sin2θcos2φ).

图2 偏振康普顿散射示意图Fig.2.Schematic diagram of polarized Compton scattering.

对出射光相互垂直的两个偏振方向x'和y'的分截面求和得到入射光为线偏振光且不测量散射光偏振的实验中的截面[25]为

一般地,如果入射光是线偏振度为PL(0

式中,pp 表示部分偏振(partially polarized),相应的图2 中X轴方向取为入射光透过理想偏振片后的透射强度最大方向.注意到公式中k/k01/[1+k0(1−cosθ)] 与 方位角φ无关,所以公式对φ积分后得到截面沿散射角θ的分布与线偏振度PL无关,进而总截面与PL无关.

下面讨论(14)式的最小值,做变换xk/k0,则

式中,x∈[1/(1+2k0),1].记(15)式右边方括号内的项为函数g(x,φ),则

∂g/∂x在x1/(1+2k0) 处的值为

令f(k0,a)0,得到k0关于a的唯一解,记为k0K(a).易知,当f(k0,a)≤0 时,方程h(x,a)0 在定义域内存在唯一解,记为xX(a).讨论函数(15)的最小值得到如下结论:若k0

由于K(1+PL)≥K(1)≈1.74,荧光CT 所用X 射线能量一般小于0.2 个电子静能(100 keV),远小于K(1),所以最终得到截面最小值点对应的θmin:

不同入射光能量下的截面最小点θmin随偏振度P的变化如图3 所示,可见荧光CT 所用的X光能量范围内,微分截面最小值对应散射角在90°附近.偏振度越高且入射能量越低,则θmin越接近90°.由于部分偏振光可以统计描述为完全偏振光(PL1 )和一个完全非偏振光(PL0)的叠加,所以最小微分截面方向接近入射光完全偏振成分的电矢量方向.

图3 不同入射能量时,θmin 随线偏 振度变化Fig.3.θmin varying with polarization for different incident energy.

图4 描述了不同入射能量下,微分截面最小值随偏振度的变化.可见在荧光CT 成像的能量范围内,随入射光偏振度的增大,最小微分截面值减小.这为荧光CT 去散射提供了思路,即入射光偏振度越高,在θθmin,φ0 or π 处的康普顿散射本底越低.

图4 不同入射能量时最小微分截面随偏振度的变化Fig.4.Minimum differential cross section varying with polarization for different incident energy.

前面讨论的情形是单个自由电子与入射光的康普顿散射截面,这种情形很理想,方便分析.实际上物质中的电子都不是自由的,此时,一个原子的总散射截面与其各个电子在自由状态下散射截面的总和一般不相等,引入非相干散射函数S(q,Z),则一个原子的非相干散射微分截面为[29]

式中,Z是原子序数;xsin(θ/2)/λ0,

λ0是入射光子波长(单位为Å).非相干散射函数的理论计算方法和数值表可在相关文献中查寻[30−33].非相干散射函数满足 0 ≤S(x,Z)/Z≤1,对于低Z原子的散射,只有在x较小(软X 射线或小角度散射,x<1 )时,S(x,Z)/Z才显著小于1.例如100 keV 入射光子,波长约为0.124 Å,入射到氧原子(Z8 ),散射角满足θ<4◦才使S(x,Z)/Z <0.5,实际实验中,入射光束总有发散角和宽度,在小散射角处测量荧光信号存在很强的入射光背景,所以实验设置应该在较大角度处测量荧光信号,因此在这个测量角度范围内单个自由电子微分截面最小值的分析依然是有效的,即考虑非相干散射函数后,入射光线偏振优势方向是依然是近似的散射强度最小点.

单个自由电子的Klein-Nishina 散射微分截面公式在kk0→0 时得到单个自由电子的Thomson散射截面公式:

若入射光线偏振度为PL,则相应的散射截面公式为

由(22)式易知,最小微分截面方向为θπ/2 且φ0,π,所以最小微分截面方向为入射光完全偏振成分的电矢量方向.从经典电动力学角度出发也可以得到相同的散射截面公式,即电子在入射光电场驱动下振动辐射电磁波.

一个原子的相干散射截面为[25]

其中F(x,Z) 是原子形状因子,可通过文献[30−32]查询原子形状因子和理论计算方法.原子形状因子随x的变化趋势与非相干散射函数相反,x越大F(x,Z) 越接近零,相反地,x趋于零时F(x,Z) 趋于Z.

根据(20)式和(23)式,一个原子的康普顿散射与相干散射微分截面的总和为

(20)式、(23)式和(24)式忽略光的偏振对原子形状因子和散射因子的影响[26],所以根据(14)式和(22)式的最小值的讨论,一个原子的总微分截面的最小值点一定满足φ0,π.

一个分子的散射截面角分布为[34]

其中m(x) 表示分子干涉函数;s(x) 和f(x) 表示采用独立原子模型(independent atomic modelling,IAM)计算得到的非相干散射函数和分子形状因子[33],即

其中Fi(x,Zi) 表示分子中第i种元素原子的形状因子;Si(x,Zi) 表示第i种元素的非相干散射函数;ni表示分子中第i种元素的原子个数;m(x) 表示分子干涉函数[34].

计算得到水分子在方位角φ0,π 处的微分散射截面随散射角θ的分布如图5 所示.从图5(a)可以看出,相干散射微分截面随散射角增大快速减小,在大角度处,微分散射截面由主要由非相干散射贡献.由图5(b)中总微分截面角分布随线偏振度的变化可知,线偏振度越高,垂直方向散射的微分截面越小.PL0.5 时不同入射光能量下的总微分截面如图5(c)所示,可见在XFCT 的光源能量范围内,不同能量光子的水分子散射截面最小值均在入射光线偏振优势方向附近,且该方向上微分截面相近.

图5 φ=0,π 处微分截面随散射角变化 (a) PL=1 时相干、非相干和总微分截面;(b) 不同线偏振度时的总微分截面;(c) PL=0.5 时不同入射能量下的总微分截面Fig.5.Differential cross section varying with θ at φ=0,π :(a) Incoherent,coherent and total differential cross section at PL=1 ;(b) total differential cross section varying with θ for different PL ;(c) total differential cross section varying with θ for different incident energy at PL=0.5.

4.2 散射光的偏振

联系康普顿散射前后光束Stokes 矢量的变换矩阵T为[25]

其中k0k0n0,kkn,S是电子初始自旋方向,n0是入射光动量方向,n是散射光动量方向.在各向同性材料中,电子初始自旋方向的平均值为0,所以矩阵T中含S的矩阵元为零.在变换矩阵T成立的坐标系中要求垂直散射平面方向完全线偏振光的归一化Stokes 矢量为即y和y'轴垂直于散射平面,由图2 易知,坐标系xOy可 由XOY绕OZ轴旋转ψ−φ得 到,则 由(4a)式、(4b)式和(4c)式可知,在XOY坐标系内Stokes 矢量为的入射光转换到坐标系xOy后归一化Stokes 矢量变为

在各向同性材料某点发生康普顿散射后,散射光Stokes 矢量变为 [If,ζ′]T,即:

由(29)式亦可得到散射光偏振的主方向与散射平面法向夹角β,符号规定为:迎着光子动量方向,在表示散射光偏振的局部坐标系x′Sy′内,y'轴以最小角度旋转到e',若其旋转方向为顺时针,则β为正,否则为负.规定β的符号后,在局部坐标x′Sy′旋转β角得到的新坐标系中,散射光Stokes矢量的第二个分量应达到正的最大值,据此条件可求得

其中函数 a tan2(y,x) 为双参数反正切函数.

在φ0 or π 散射平面上散射光的归一化Stokes矢量、线偏振度和偏振主方向角为β分别为

特别地,在θπ/2,φ0,π 散射方向上散射光的归一化Stokes 矢量和线偏振度为

图6 展示了在不同入射光能量下φ0,π 散射平面内随散射角θ的变化.图7 展示了不同能量下θπ/2,φ0,π 散射方向上散射光线偏振度随入射光线偏振度的变化.若入射光为完全线偏振光(1 ),则θπ/2,φ0,π 散射方向上散射光的线偏振度为零,这与(12)式是相符的.从图6 可以看出,在XFCT 所涉及的能量范围内,相同散射角处不同能量下散射光线偏振度差距较小,而在θπ/2 附近,部分偏振的入射光,其散射光也具有较高偏振度.从图7 可知,只有当入射光线偏振度接近完全线偏振时,θπ/2 的散射光的偏振度才显著降低.

图6 φ=0 or π 处 随散 射角变化Fig.6. varying with θ at φ=0, π.

图7 不同能量下 θ =π/2,φ=0, π 处散射光线偏振度随入射光线偏振度变化Fig.7.The linear polarization of scattering photons varying with at θ =π/2,φ=0, π for different incident energy.

从(30)式、(31)式和(32a)式可知:完全非偏振入射光的散射光的偏振主方向始终垂直于散射平面,部分偏振入射光的散射光偏振;部分偏振入射光在φ0,π 平面内的散射光偏振度是散射角的函数,且存在两个零点,左零点左侧和右零点的右侧偏振优势方向平行于散射平面,在两个零点之间,偏振优势方向垂直于散射平面,在两个零点处散射光完全非偏振;在φ0,π 散射平面内,完全线偏振入射光的散射光在θπ/2 时的偏振优势方向均平行于散射平面,在θπ/2 时散射光完全非偏振.

根据以上关于散射光的线偏振度及荧光线偏振度的分析,并假设存在理想的X 光偏振片,因为K 层荧光的线偏振度为零,K 层荧光经过理想偏振片后减少50%,而散射光具有一定线偏振度,比如当入射光完全非偏振时,在单个自由电子这个理想的情况下,垂直于入射光方向的散射光的线偏振度可达到90%,经过理想偏振片后散射光可减少95%,因此,在经过偏振片后,虽然荧光强度减半,但是荧光强度与散射强度比值即信噪比在理想情况下可提高十倍,所以在原理上,散射光的线偏振度和K 层荧光线偏振度的差异可用于提高荧光CT 成像的信噪比.

5 结论

本文分析了XFCT 的物理过程及其性质.计算荧光退激发时间和荧光产额发现,中高Z元素荧光产额基本饱和,而荧光退激发时间远远小于现有探测器时间分辨率,因此,荧光退激发时间无法用于分辨散射光与荧光,提高原子序数可以增加荧光的能量,荧光能量越高,其穿透性越强,有助于提高探测到的荧光光子数,减小统计噪声,提高信噪比.偏振光源有助于减少康普顿散射背景,入射光偏振度不影响康普顿散射总截面,但对截面的角分布有影响,部分偏振入射光在偏振方向附近具有最小的散射截面,在此方向上探测荧光光子有助于降低康普顿散射背景,并且对于完全非偏振或部分偏振的入射光,在大部分方向上的散射光仍然具有一定的线偏振度,而K 层荧光线偏振度为零,因此荧光和散射光线偏振度差异可以用于进一步降低荧光探测中的康普顿散射背景,由于实际成像过程的复杂性,通过偏振性质来减少散射背景需要进一步模拟和实验验证.