张智河 陈顺云 刘培洵 刘琼颖
1)首都医科大学,生物医学工程学院,北京 100069 2)中国地震局地质研究所,地震动力学国家重点实验室,新疆帕米尔陆内俯冲国家野外科学观测研究站,北京 100029
近年来的基岩温度观测证实: 通过观测基岩温度可以获取地壳应力的动态变化信息,即所谓的 “热测应力”(Chenetal.,2016,2019;陈顺云等,2020)。陈顺云等(2020)将基于康定MS6.3地震同震温度响应获取的同震应力变化量级和空间分布特征与测震学方法得到的结果进行对比,发现两者一致,证实了根据野外温度观测分析同震应力变化的有效性。将来,伴随着温度测量技术的进一步提高,热测应力有望为地震研究带来新的机遇。
上述野外观测(Chenetal.,2019;陈顺云等,2020)所采用设备的温度测量精度为0.2mK(张智河等,2018),温度分辨率为0.01mK。另外,杨晓明等(2018)采用将20个Ptl000串联的思路,基于惠斯通平衡电桥设计了高分辨率温度传感器。经测试,其温度分辨率的量级可达0.01mK,测量精度不详。根据文献数据(杨晓明等,2018)可知,上述仪器的测温精度应优于1mK量级,考虑到其测量的是室内静置水桶内水温,难以满足高精度温度标定的要求,故仪器实际的测温精度可能更优。总体而言,目前温度的观测精度依然有限,限制了热测应力的深入应用。
特别地,地震引起的应力变化: 1)与测点至震源的距离有关,且随距离的增加快速衰减;2)与震级有关,随震级下降呈指数衰减,如在上述康定MS6.3地震中观测到的同震响应在随后的MS5.8地震中并没有出现(陈顺云等,2020)。欲使热测应力在地震研究中发挥更大的作用,尤其是探测强震前的应力变化信息,亟需研发更高精度的温度测量系统。
本文拟在原温度观测系统的基础上,从传感器和传输线效应以及抑制系统噪声等方面改进温度测量精度,以期从观测技术上推动热测应力相关研究向前发展。
原有的观测精度达0.2mK的野外观测设备(张智河等,2018)采用的是恒流源激励的两线制方式,主要存在以下几个方面的问题: 1)恒流源激励的两线制,其温度信号采集点在电路板上,测温回路包含传感器、引出电缆及模拟开关的电阻,同时长导线使测温远端相对电路板端存在温差热电势;由低温漂定值电阻提供的负差分输入端在电路板上,温度测量可视为对地的单端测量;2)单端温度信号经模拟开关直接引入模数转换(AD)芯片,再由AD芯片内部的可编程放大器放大。低温漂定值电阻提供的参考电压相对基准电压的一半偏低,不能充分利用AD芯片的满幅测量范围;3)温度计算模型过于简单,带宽抑制能力有限,采集处于铂电阻自加热升温段。
为了提高测量精度,从以下几个方面进行了完善: 1)优化传感器及其传输线;2)抑制电子电路噪声;3)消减温差引起的热电势效应。以上3个方面均会影响温度的测量精度。采集系统的基本原理如图1 所示,包括几个关键步骤: 1)换向恒流桥式测温;2)前置放大模拟滤波;3)应用现代数字滤波技术进行数据前处理,尽可能使得数据采集处于铂电阻自加热(准)平衡时段。下文拟分别对其进行考虑。
图1 恒流桥式温度测量功能框图Fig.1 Functional sketch map for bridge-type temperature measurement system based on the constant-current source.
两线制测温的驱动线与信号线混用(路傲轩等,2018),测温分支点在电路板模拟开关处。因导线较长,有的通道长达几m,使得导线电阻的温度效应、电流流过传输线的电压效应及长导线温差电势效应都混合在信号中,且属于差分干扰信号,均混合于所采集的数据中,因此不能抑制温度传感器之外的混杂信号。
为了突破两线制的局限性,本文引入四线制桥式测温,采用惠斯通平衡电桥,将定值参考电阻与测温电阻对称放在测温端。通过野外观测数据可知温度的波动幅度<0.1K,可视为恒温环境。温度的灵敏度受传感器金属外壳热容量的影响,其热容量较大有利于滤除温度信号的高频成分。由于电路板上的温度环境复杂,将定值参考电阻置于测量端较为合理。将四线制测温驱动线与信号线分开,测温分支点设置于温度传感器处,信号输入高阻抗,信号通路的电阻变化可忽略。桥式测温取自测温电阻与定值电阻的差分信号,信号地干扰信号相对差分温度信号输入为共模信号,可由差分放大系统的高共模抑制比抑制。恒流与恒压相比其供电温度信号的线性度更好,适合较宽的温度范围测量,是高精度测温的首选。
当连接不同材质金属导体的连接处存在局部温差时,其两端就会产生一个热电动势(热电势),该电动势的方向和大小与导体的材料及接点的温度有关,该现象被称为热电效应。热电势在测温回路中是普遍存在的,其温度系数与测温电阻不同,且非线性特征明显。
为了消减热电势对传感器的测温影响,利用在同一温度下热电势的方向和大小不变的原理(丁炯等,2018),采用正、反方向加电消除。设温差热电势为VT,桥臂电流为I,测温电阻Rs=Rs1=Rs2,定值电阻Rd=Rd1=Rd2,如图2 所示。当电流为正向时,
V1=I·(Rs1-Rd2)+VT
(1)
当电流为反向时,
V2=I·(Rs2-Rd1)-VT
(2)
取二者平均值,即式(1)+式(2),有:
(3)
式(3)中已不包含温差热电势VT。
图2 恒流桥式温度测量电路Fig.2 Temperature measuring circuit on the base of constant current source and bridge-type sensors.
理论上,与两线制相比,四线制有诸多理论上的优势。惠斯通电桥测温是传统的测温方法之一,恒流源驱动能满足高线性、高精度的要求。然而,如后文所示,相对两线制而言其改善不算明显,可能的原因是模拟放大倍数不够,通过信号不足以分辨出不同制式的区别,采集得到的数据的噪声水平均接近系统本底噪声,制式变化信息完全淹没于其中。另外,恒流源正、反向2次测量可有效抑制长导线温差的热电势效应,同时采用斩波稳零运放可有效抑制电路低频1/f噪声的作用,属于模拟信号处理的范畴。
实际上,将模拟信号转换为数字信号后,还可以采用数字处理技术对噪声进行抑制。本文拟引入Kalman滤波法对白噪声进行压制。本文中关于噪声抑制处理可简单地分为2部分,即模拟信号降噪和数字信号处理。
1.3.1 模拟电路降噪
恒流桥式温度信号的输出与激励电流呈线性关系:
Vout=A·I·(Rs-Rd)=A·I·R·α
(4)
其中,A为放大倍数,I为分支电流,α为温度系数。
电阻热噪声为
(5)
其中,k为玻尔兹曼常数,T为绝对温度,BW为带宽。
信噪比为
(6)
噪声与带宽相关,前置差分放大器、低通滤波器可有效限制带宽,增强AD的抗混叠能力,AD数据经过内部sinc3滤波后输出,低输出速率能保证高有效位数的转换精度。高的基准参考电压也能相对增加信噪比。
1.3.2 数字信号滤波
将模拟信号转换为数字信号后,采用Kalman滤波方法进一步压制白噪声。Kalman滤波法于1960年提出,是基于维纳滤波法的改进方法,通过递归方法克服了维纳滤波褶积运算的缺点。同时,Kalman滤波法也是贝叶斯预测的一种特例,故其预测具有无偏、稳定和概率意义上最优的特点,尤其适用于线性系统中以抑制高斯噪声(杨超等,2020)。
本文中的电阻性测温系统相当于一维情况,所采集的数据c包含测量量x及测量噪声v,v的均方差为R,测量方程为
c=x+v
(7)
其中,x为状态量,w为状态噪声,w的协方差为Q,无控制量的状态方程为
y=x+w
(8)
预估值协方差为P,当前估计值为x,采集值为c,根据Kalman滤波公式推算新的预估值y:
y=x+kg(c-x)
(9)
其中,Kalman增益为
(10)
新的预估协方差为
P=(1-kg)×(P+Q)
(11)
将式(11)代入式(10)以计算新的Kalman增益kg。设定Q、R,预设1个偏离的初值P,经反复迭代得到P和kg序列。由式(9)可以看出,Kalman滤波是一个滑动比例压缩的过程,对偏差大的数据压缩度更大,由此使数据列的离散度变小,提高了精度,从滤波效果看属于低通滤波,但存在相位偏移。对于时变信号而言,其幅值将被压缩,对于准温度平衡信号而言则可认为其不受影响。
由式(10)可知,kg最终趋于稳定的<1的值,Q、R、P的初值由人为设定,根据实测数据优化选择其数值。理论测量精度的改善程度与kg相关,Kalman滤波是否收敛及其收敛程度受实际观测数据和采集时长的影响。由于单片机的运算能力与低功耗约束的限制,未采用自动参数设定功能。
温度测量采用测温电阻和定值电阻组成对称桥的方式实现,如图2 中的虚线框部分所示。其中,Rs1、Rs2为测温电阻,Rd1、Rd2为定值电阻,2个桥路均由测温电阻与定值电阻串联组成,其阻值相同,通过的电流为恒流源总电流的一半,对称桥式结构不因电流源换向而改变支路的电流分配。
在选择温度测量电阻时,设计平衡温度点为15℃,定值电阻为2.117kΩ,温度系数为1ppm,测温电阻由2只Pt1000串联组成。桥臂电流为0.8mA,定值电阻提供的参考电压为1.694V,略高于AD参考电压源(3V)的一半,抬高了单电源前置放大平衡点的电压(图3),信号经放大后即可充分利用模数转换范围。设计的温度范围为-15~45℃,前置放大器的增益为16倍。AD采用24位Σ-Δ温度采集专用芯片AD7714AR-5,自带放大功能,本文中设置增益为1倍,温度分辨率为3.58μK。
图3 温度测量的前置放大器与模数转换电路Fig.3 Circuits for preamplifier and analog-to-digital conversion of temperature measurement.
温度测量模块有4组测量通道。如图2 所示,K1、K3分别用于切换恒流驱动源与测量通路,K2用于切换恒流源驱动方向。定值电阻起参考电压的作用,Rs1和Rs2分别反映换向前、后的温度变化。TC1用于降低电流源通路的噪声,TC2用于降低信号输出通路的噪声,起到低通滤波的作用。
前置运算放大器隔离缓冲输入信号,阻容低通滤波电路限制模拟信号带宽,满足AD转换的带宽约束条件,电路如图3 所示。
温度采集配置的转换速率为50Hz。由AD触发采集中断,采集缓冲105组数据,在采集数据的同时进行Kalman滤波。
设置初值Q=0.04、R=1.0、P=1;利用式(10)、(11)迭代计算Kalman增益kg。定义浮点数组,其有效位数与24位AD的位数相同,采用32位自定义浮点格式: D31为符号位;D30为0;D29—D24为指数位,负数用补码;D23—D0为小数位,均为正数,首位均为1。限于24位精度取43组数据保存在固件Flash中。
读取所采集的第1组数据作为基准数,并用后续采集的数据减去基准数后转化为自定义浮点数。从第2组数据起,将后续数据代入Kalman滤波式(9)中,设x的初值为0,查表读取第1组Kalman增益kg,计算新的最佳估值并更新x,进行第1次Kalman滤波。将此最佳估值代入Kalman滤波公式进行第2次Kalman滤波,第2次滤波时x1的初值为0,kg不变,计算新的最佳估值并更新x1,依次读取采集数据与Kalman增益kg计算新的最佳估值。第3次Kalman滤波的效果改善不大,故仅取到第2次滤波的结果。经43次迭代后,将得到的第2次最佳估值转化为整数作为修正数,与基准数代数求和,更新为新的基准数。修正数再次转化为浮点数以修正第1次、第2次最佳估值,以最后1组Kalman增益kg反向迭代到起始采集数据处,进行定值Kalman滤波。
再从采集数据起始处进行2次正向定值Kalman滤波,对第2次最佳估值进行两两求和,得到50组数据。
电流换向,信号也被反向,对采集数据进行补码运算,保持与正向采集数据同相。进行与上述过程相同的处理,并进行2次Kalman滤波,对第2次最佳估值两两求和,并与正向滤波得到的最佳估值求和,得到50组数据。
利用与上述过程相同的方法将50组数据按照逆向顺序进行2次Kalman滤波,修正相位偏移,对迭代的第2次最佳估值求和后除以200;正、反向基准数取平均值后与滤波得到的最佳平均估值进行代数求和得到最终测量值。通过同样的方法循环采集得到4个通道的最终测量值。
由于Q、R值的设定受人为因素影响,迭代算出的最终Kalman增益kg是滤波效果的决定参数,滤波效果同时受具体采集数据的影响,下面给出一定值Kalman滤波实际算例以供参考,其中kg=0.176i556。
图4 给出了温度采集的Kalman滤波效果。正向采集的均方根误差为40.05,第1次Kalman滤波后的均方根误差为29.93,相对原始数据改善了1.34倍;第2次Kalman滤波后的均方根误差为27.13,相对原始数据改善了1.48倍。反向采集的均方根误差为11.72,第1次Kalman滤波后的均方根误差为5.87,相对原始数据改善了2.00倍;第2次Kalman滤波后的均方根误差为4.43,相对原始数据改善了2.64倍。对上述正、反向采集的第2次Kalman滤波后的数据取均值,其均方根误差为13.94,再对均值数据进行2次Kalman滤波,第1次Kalman滤波后均方根误差为12.52,第2次Kalman滤波后均方根误差为10.28。综合考虑原始数据的均方根误差为29.51可知,最终结果相对改善了2.87倍。
图4 温度采集Kalman滤波的效果Fig.4 Results of Kalman filtering to temperature acquisition.a 正向采集;b 反向采集;c 正、反向采集(滤波后)的结果
在室内实现温度波动<0.01mK的测试环境极具挑战,需要投入大量的人力物力。前期经观测发现基岩的温度极为稳定,本文采用野外观测数据进行仪器评测,以评估仪器的测量精度。同时,还可以验证仪器对于野外工作的适用性。具体结果如图5 所示。
图5 不同版本仪器的基岩温度观测结果Fig.5 Results of bedrock temperature observation from different versions of instrument.a 新疆测点结果(2019年,四线制Ⅱ);b 云南测点结果(2017年,四线制Ⅰ);c 云南测点结果(2017年,两线制)。所有温度测点均已去除趋势变化
从图5b、c可以看出,两线制的温度测量精度达0.2mK,四线制的温度测量精度达0.1mK,峰-峰值分别为-0.44mK和0.47mK。与两线制相比,四线制的温度测量精度提高了1倍。
四线制比两线制具有诸多理论上的优势,惠斯通电桥测温可充分提取温度差分信号,本研究的野外观测结果图5b所示,达到0.1mK量级,但四线制与两线制相比所具有的优势并没有明显地体现出来,推测可能是由于所采集的数据噪声水平接近采集系统本底噪声,放大倍数不够,信号分辨率不足以反映出制式不同的区别。
在上述四线制的基础上,重点利用现代数字滤波的处理技术抑制采集信号的噪声,以提高信噪比。同时,改进电源管理以不断提高模拟电源纹波抑制能力,运用滤波电路限制信号带宽、抑制AD混叠。理论上经过这些抑制措施后,可以将测量精度提高3.76倍。
图5a给出了改善后的四线制采集板的实际观测结果。从图中可以看出,温度测量精度达0.03mK,峰-峰值分别为-0.13mK和0.12mK。相比未经电路板噪声抑制的四线制采集板测量结果(图5b),新方案的温度测量精度提高了3倍,与理论预期比较接近。
本文基于低温漂定值电阻与测温电阻组成的平衡桥式四线制温度传感器,考虑了恒流换向驱动和Kalman数字滤波等一系列技术改进措施后,成功研发了新一版高精度温度测量系统。经过野外观测检试,取得以下成果:
(1)最新一版温度观测系统在温度测量精度上得到明显的提升,温度分辨率为0.003mK。经野外观测数据测试,其精度达0.03mK。
(2)Kalman滤波方法对于抑制电子电路的噪声水平有明显的优势,本文中将精度提高了约3倍,可为高精度温度测量仪器的设计提供参考。
(3)未经噪声抑制的四线制测温系统,其温度测量精度达0.1mK。与二线制测量系统相比,其温度测量精度提高了1倍。
库仑应力变化对于地震的触发作用一直是地震研究的热点之一,如何获取同震库仑应力变化值一直是该领域的难点问题,热测应力可为观测库仑应力变化提供一种新的选择。同震库仑应力变化的量级为0.01~0.1MPa,由于地壳典型岩性的热应力系数约为1mK/MPa(Yangetal.,2017;Chenetal.,2019),观测上述量级的应力变化,要求温度测量精度达到0.01~0.1mK量级,这是一项极具挑战性的工作。
综上所述,最新一版的温度测量系统的整体性能已获较大幅度的提升,测量精度已达0.03mK,从技术上可以获得0.03MPa的动态应力变化,达到了同震库仑应力变化测量的量级。野外观测证实了其测温方案具有可行性,从测温技术层面较大幅度提高了热测应力的可测量范围,这对于推动热测应力研究的发展大有裨益。
致谢审稿专家为本文提供了细致的修改建议,在此表示衷心感谢!