燕芳，毛莉程，张俊林

(内蒙古科技大学 信息工程学院， 内蒙古 包头 014010)

由太赫兹波穿透样品后的时域脉冲信号可得到太赫兹波磁场的幅频与相频信号，太赫兹时域光谱系统(THz-TDS)可以探测出此脉冲信号并进一步处理得到该样品的折射率谱和吸收谱，即样品的色散及吸收等特性信息，因此可以利用太赫兹技术进行物质的定性及定量分析[1].

由于固态样品分析中样品颗粒尺寸与太赫兹波波长接近，故会造成太赫兹波透射中的强烈散射问题.由THz-TDS获取的吸收谱包括2部分，即样品自身对太赫兹波的本征吸收谱αabs，和由散射带来太赫兹波衰减的散射吸收谱αsca.目前，定量分析方法直接将总吸收谱(αabs与αsca之和)作为样品的吸收谱，而忽略了散射吸收谱αsca.因此，随着太赫兹定性、定量分析准确性要求的提高，去除或降低散射效应的影响，获得样品真正的吸收谱就显得十分必要[2].

随着太赫兹技术的发展，以及太赫兹谱的应用日渐增多，对散射效应地研究已成为太赫兹光谱领域的一个重要课题.德国Freiburg大学的FRANZ M等认为固体颗粒物质的太赫兹吸收谱中存在的基线与样片颗粒对太赫兹波的散射吸收有关.澳大利亚阿得雷德大学的ABBOTT D教授团队在系列文献中利用多种计算方法及理论模型针对不同情况下的样品在太赫兹波段的散射效应进行了研究，值得注意的是，他们的方法在对太赫兹波不吸收的致密度较高的样片实验中取得了显著的效果.纽约新泽西理工大学的BANDYOPADHYAY A等人利用Mie氏散射理论研究了在松散分布情况下的硝酸铵、面粉以及盐粒在太赫兹波段的散射效应.英国利物浦大学的SHEN Y C等人研究了颗粒尺寸在50～250 μm范围内的对太赫兹具有吸收特性的样品和不吸收样品，给出了对太赫兹不吸收物质(如PE)的散射效应与频率、颗粒度的经验表达式.哈尔滨工业大学的刘文军等人对存在多重散射的介质进行了研究，他们认为多重散射场下太赫兹时域脉冲符合不同的延迟时间相关的统计分布[3-5].

综上所述，国内外学者利用多种算法及理论模型研究了太赫兹波的散射效应，但其中大多数方法都需要预先掌握样品的参数信息，如样品厚度，各组分的折射率，颗粒的尺寸和形状等，不利于对散射效应的在线快速测试计算，而且多采用对太赫兹吸收不明显的样品，如PE等.本文的基线拟合模型以样片从THz-TDS系统获取的实验数据为基础，适用于对太赫兹波有尖锐吸收的样品(氨基酸样片)，且可适用于样片厚度及各组分折射率等样片参数未知的情况，实用性较强.

1 散射基线拟合模型的建立

1.1 吸收谱基线的产生

以氨基酸的太赫兹时域光谱测试为例说明基线的产生原理.实验中采用透射式太赫兹时域光谱仪(THz-TDS，EKASPLA公司)，飞秒激光器是钛蓝宝石锁模激光器(Spectra-Physics公司)，它产生的飞秒激光照射到砷化镓(LT-GaAs)晶体上可以产生太赫兹脉冲.整个太赫兹光路系统在实验中会始终处于充入氮气的密闭试验箱内，箱内温湿度均经过严格控制以减少空气中水分对太赫兹波的吸收.

实验用纯氨基酸样品(谷氨酰胺、组氨酸与苏氨酸)均购于WakoPure Chemicals公司，由于纯样品难以成型，故将其与聚乙烯粉末(PE,其在太赫兹波段几乎透明)混合均匀，在玛瑙研钵中研磨后用10 MPa压成直径约13 mm，厚约1.0 mm的圆片.将样片置于透射型太赫兹时域光谱系统中获取该氨基酸样品的吸收谱.

氨基酸样品在太赫兹波段的介电参数提取方法具体参见Duvillaret[6-7]等文献.将透过样片的太赫兹信号记为Esam(ω)，没有透过样片在空气中传播的太赫兹信号称为参考信号，记为Eref(ω).这2个信号是由THz-TDS系统测出的时域内的样品信号和参考信号取平均值后再进行快速傅里叶变换后得到的，由二者可以计算出样品的吸光度α(ω)，如式(1)：

(1)

式中：D为样片厚度，mm；n(ω)为样品的折射率.

图1为谷氨酰胺样品的太赫兹吸收谱，图中“基线”(虚线所示)是由于样品颗粒对太赫兹波散射造成的，基线的存在改变了吸收峰的中心频率和峰宽，会在后续对样品的定量分析中引入误差.

样品颗粒大小是影响散射基线斜率的主要因素，故利用式(1)将不同粒度的组氨酸和苏氨酸样片在太赫兹波段的吸收系数计算出来，如图2所示.

对比图2中2种氨基酸样片在粒度不同情况下的吸收谱可见，二者的吸收谱线型相似，峰型相似，吸收峰所对应的频率值也相同，区别在于峰高不同.粒度大(<150 μm)的氨基酸样片吸收峰的峰高大于粒度小(<57 μm)的氨基酸样片.可见，当样品颗粒大小接近入射波的波长时，样品的散射效应更为明显.此时，入射的太赫兹波穿过样品时损失的部分应包括样品本身对太赫兹波的吸收，以及由于样品散射吸收而造成的衰减，故会在粒度大的氨基酸样片吸光度曲线上造成斜率较大的“基线”，它会使吸收谱整体产生畸变，抬高吸收峰，甚至改变吸收峰的中心频率和线宽.实验结果说明散射光强会随着样品颗粒粒度增大而增强，吸收谱中的基线斜率也随之增大.

1.2 散射基线拟合算法原理

由图2实验结果对比分析可知，样品颗粒的粒度越大，散射光强越强，吸收谱中的基线斜率越大.MASON W P等人也已证实当颗粒尺寸增大到可以和入射光的波长相比拟的时候，散射效应近似的同频率的平方成正比.基于以上研究基础，由于样品的本征吸收与散射吸收与频率具有不同依赖关系，结合样品粒度构造出一个基于数据及粒度因素的表征散射吸收的模型，利用该模型可将吸收谱中的本征吸收峰信息从散射吸收背景中提取分离出来，形成表征物质结构特性的本征吸收谱，利用该模型拟合基线可以减小吸收谱中的散射损耗带来的定量分析误差[8].

1.2.1基线拟合建模原理

当一束太赫兹波信号垂直入射到复介电常数为n(ω)+ik(ω)，厚度为D的样片上后，把透过样片的THz信号的电场强度经傅里叶变换后记做样品信号Esam(ω)，没有透过样品而是在真空中传播距离L后的场强信号做傅里叶变换后记做参考信号Eref(ω).

(2)

式中：α(ω)为样片总的吸收系数.由朗伯定律吸收系数α(ω)与消光系数即复介电常数的虚部k(ω)有如下关系式：

(3)

由式(2)可求出对应的振幅谱如下：

(4)

如前所述，由THz-TDS获取的吸收系数α(ω)为样品自身对太赫兹波的本征吸收系数αabs(ω)与由样片散射引起太赫兹波衰减的散射吸收αsca(ω)之和，即：

α(ω)=αabs(ω)+αsca(ω) .

(5)

(6)

由实验结果分析可知，样品颗粒粒度与样品散射吸收成正比.故而，构造出一个表征散射吸收的模型，将散射吸收表征为角频率ω的指数函数，在指数中引入样品的颗粒直径d与波长的比，如式(7)所示.

(7)

式中：α0为样片本身参数相关的造成散射效应的待定系数，可以称为散射系数，cm-1；d为样片颗粒直径，mm；λ为入射太赫兹波的波长，m；ω0为αsca=α0时对应的角频率，rad·s-1.将αsca带入式(6)，得出样品信号的振幅谱如下：

(8)

(9)

由式(9)可知，散射系数α0不会改变Z(ω)的谱线形状，只是使其谱线上下平移.为了消除这部分由散射引起的谱线平移，我们令指数A中的系数A从0开始以0.01的步长增加，迭代计算式(10)中的Z0(ω)，其中Z1(ω)为Z(ω)的均值，迭代计算直至Merr的值达到最小为止，最后，根据Z0(ω)的值就可以计算出表征样品本征吸收特性的吸收系数，如式(12).将由THz-TDS获取的吸收系数与相减即可得到基线.

Z0(ω)=Z(ω)-Z1(ω) .

(10)

(11)

(12)

1.2.2利用基线拟合修正吸收谱

利用基于数据处理并结合粒度因素的散射模型拟合基线，修正吸收谱的具体步骤如下：

(1)由THz-TDS系统获取同一环境和温度下的样片参考时域谱Eref(t)及信号时域谱Esam(t)，经傅里叶变换得到Eref(ω)及Esam(ω)，进而计算出原始吸收谱α(ω)；

(2)利用基于可变阈值的峰值探测算法结合测量系统的频率分辨率确定拟合基线频率区域ωS；

(3)令系数A从0开始以0.01的步长增加，迭代计算式(9)～(11)，直至Merr的值达到最小；

(4)将步骤(3)中得到的A值代入式(12)，计算出本征吸收系数αabs(ω)；

(5)由式(5)计算出散射吸收并在吸收谱曲线上绘制出基线.

1.3 基线模型中关键参数的确定

在上述拟合算法中，拟合基线的频率区域ωS的选取是难点问题.在结合实验装置中频率分辨率的基础上，采用峰值探测算法来解决这个问题，下面以组氨酸样品为例说明解决方法.

图3为组氨酸(His)的吸收谱.首先，在获取样品吸收谱的基础上，利用峰值检测算法确定吸收谱中特征吸收峰的峰位(图3中圆点标示位置)，其次，根据测量计算吸收谱时的频率分辨率来确定提取基线频率区域ωS(图3中交叉号标示位置)，该频段处的吸收谱线可以认为是单纯由样片散射引起的太赫兹波吸收损耗，通过在此频域内的迭代计算最终确定指数A的值.上述迭代计算的步长原则上可以任意选择，步长越小则估算精度越高，但步长过小会加大运算量，降低运算速度，故应综合考虑选择合适的迭代步长.

2 基线拟合模型仿真结果分析

选取组氨酸(L-Histidine, His)和苏氨酸(L-Threonine, Thr)作为测试样品(质量分数为99.0%).选用购于Sigma-Aldrich公司的聚乙烯(polyethylene, PE)作为样品稀释剂，PE晶体粉末的颗粒尺寸约为50～75 μm.实验前将氨基酸粉末过标准筛，将其粒度控制在小于57 μm(260目筛)与小于150 μm(100目筛)两种情况下.把过筛后的氨基酸颗粒与一定量PE都加入玛瑙钵内进行研磨混合，待混合均匀后用8 MPa的压力压成约1 mm厚，约13 mm直径的圆形样片.实验中按上述方法制作了不同粒度与折射率的组氨酸与苏氨酸样片共7片，其中6，7号样片为2种氨基酸的混合样片，样品配比见表1.

表1 氨基酸样片配比表

利用上述基线模型对表1中不同粒度氨基酸样片的吸收谱进行修正.为了比较粒度因素对散射基线拟合的影响，将同种类氨基酸在2种不同粒度(150 μm与57 μm)下的样片拟合结果绘制在一张图中进行对比.图4(a)，(b)分别为组氨酸样片与苏氨酸样片的拟合结果对比.

由图4可见，利用上述结合粒度因素的散射基线拟合算法可以拟合出吸收谱中的散射基线(如各个子图中虚线所示)，去除基线后的修正谱中吸收峰的中心频率基本不变，峰高略有减小.对于同种类的氨基酸样品来说，粒度大的样片基线斜率明显大于粒度小的样片，由此可见，粒度越大的样片，其散射损耗也越大，这点之前实验得出的结论完全吻合.

3 基线拟合模型的有效性验证

为了进一步说明这种基线拟合算法的有效性，对表1中的混合氨基酸样片(6，7号样片)采用基于最小二乘法的定量分析，通过定量分析精度说明该算法的有效性.

首先，利用传统吸收谱对表1中的混合样片进行定量分析，由传统吸收谱计算出的混合样片中各个单质氨基酸样片的摩尔浓度如表2第五列所示.然后，利用该算法将单质氨基酸样片的传统吸收谱进行修正，得到单质氨基酸样品修正后的本征吸收谱.最后，利用修正后的本征吸收谱计算出混合样片中各个氨基酸组分的摩尔浓度，计算结果如表2第八列所示.表中第四列为混合样片中的苏氨酸与组氨酸单质成分的真实浓度.通过对比基线拟合前后基于传统吸收谱与修正吸收谱计算所得的混合氨基酸样片中的各个单质氨基酸样片的摩尔浓度，基于传统谱计算的摩尔浓度相对误差如表第六列所示，基于修正谱计算的摩尔浓度相对误差如表第九列所示，并比较了二者计算下的全局误差，分别如表第七列与第十列所示.基线拟合模型有效性验证结果如表2.

其中，全局误差(global error)定义如下：

(13)

式中：Pical为氨基酸计算摩尔浓度，mol·L-1；Pireal为氨基酸实际摩尔浓度，mol·L-1.

表2 基线拟合模型有效性验证结果

如图5所示，基于提出的散射基线拟合模型修正后的吸收谱定量分析所得的各氨基酸浓度均较基于传统吸收谱定量计算所得浓度值更接近其真实值，定量分析的相对误差以及全局误差均有减小，表明用本文提出的散射基线拟合算法修正后的氨基酸太赫兹吸收谱更为精确.

4 结论

本文结合现有基线拟合思路及散射因素分析结论，提出一种基于数据及粒度因素的散射基线拟合模型.利用该模型可在样片参数部分未知的情况下，实现将吸收谱中的本征吸收峰信息从散射吸收背景中提取分离，最终得到表征物质结构特性的修正吸收谱.利用该模型拟合的基线斜率与样片粒度大小成正比，能够部分消除粒度因素对太赫兹吸收谱中散射损耗的影响.

利用该法对测得的吸收谱进行基线拟合后，由修正吸收谱进行定量分析所得的各氨基酸浓度均较基于传统吸收谱定量计算所得浓度值更接近其真实值，定量分析的相对误差以及全局误差均有所减小，进一步表明本文提出的散射基线拟合模型是有效的，基于该模型修正后的吸收谱定量分析精度均有所提高.