亚地震频段通常指的是1～6 h周期的振动信号,通过计算不同仪器的亚地震频段背景噪声水平,计算结果可以作为比较不同仪器性能的依据(邢乐林等,2010)。为了消除仪器漂移和残余潮汐信号的影响,更清晰地识别和分析重力数据中的信号、更准确地研究地球的自由振荡和其他相关现象。

在过去的研究中,使用超导重力仪进行观测的台站噪声水平可能非常低,大约为0.07 μGal。在这种情况下,如果采样间隔为1 min,并且希望探测到1 nGal的信号,那么至少需要2 000个采样点,数据长度至少为14 d。因此,选择连续15 d的最平稳时期进行亚地震频段的观测,可以确保信号的检测和分析具有高准确性和可靠性。这个时间窗口足够长,可以提供足够的数据点,从而确保信号检测不受随机噪声影响,也能够满足频谱分析需求。

本文采用与超导重力仪数据处理一致的方法,计算弹簧重力仪的亚地震频段的背景噪声水平,具体步骤如下:①选取重力台网以天为单位的重力和气压数据,不做任何处理,数据采样率为分钟采样。②扣除固体潮影响。③修正气压数据中的间断、阶跃和尖峰等问题,并进行气压改正(大气导纳值为-0.3 Gal/mbar)。④对数据使用高通滤波器滤波,截止频率设定为8 h。⑤经过上述处理后,计算每一天的重力残差,并滚动计算连续15 d的RMS值,滑动窗口为一天。⑥将RMS值由小到大排序,选取15 d的RMS值最小且连续的原始重力数据,进行傅立叶变换,得到这15 d的平均振幅谱,计算其平均功率谱密度(PSD)。⑦计算SSNM值。

潮汐频段的范围是信号周期大于6 h。通过对台站相对重力仪连续一年的观测数据进行分析,计算潮汐频段内的背景噪声水平,并且与其它重力仪的背景噪声水平进行比较,发现潮汐改正和气压改正可以显著改善噪声水平(韩宇飞等,2015)。在评估潮汐频段的背景噪声水平时,首先选取台网以天为单位的原始观测数据,数据采样率为分钟。修正重力和气压数据中的间断、突跳和阶跃等问题,将1 min采样率降为1 h采样率,利用Eterna或Venedikov方法进行调和分析(Wenzel,1996; Venedikov et al,2003)。在Eterna方法的应用中,采用Tamura引潮位展开表,并采用高通滤波消除仪器漂移,当只考虑气压和极潮对台站的影响时,可以将某一特定时刻t的观测值y(t)表示为(Tamura,1987):

式中:M代表波群的总数; 对于第m个波群中的第n个潮波分量,A_mn、ω_mn和φ_mn分别代表理论振幅、角频率和初始相位; δ_m和Δφ_m是第m个波群中的振幅因子和相位滞后,需要通过计算得到; Dr(t)表示仪器的零点漂移; P(t)和Pole(t)分别是气压和极潮的观测值; a_k和b是气压和极潮的回归系数; ε(t)代表观测误差。

Venedikov等(2003)提出了一种潮汐分析方法,该方法利用48阶奇偶滤波器计算潮汐的周日、半日和1/3日波分量,这些分量被用来估计潮汐参数。目前,国家重力台网中心运用这种方法,结合Nakai检验方法和数字滤波技术来评估连续重力观测资料。Venedikov等(2003)引入Tamura的1 200个分波潮波表,并将气压等辅助观测序列、分段多项式方法纳入潮汐分析数学模型,增强了程序分析潮汐观测数据的能力(Venedikov et al,2003)。在Venedikov方法的应用中,可以将固体潮某特定时刻t观测值y₁(t)表示为:

式中:N为波群总数; h_j是第j个潮波的观测振幅; ω_j是第j个潮波的角频率; 当t=0时,φ_j为该时刻的初相位; P_k(t)是逼近零点漂移所选择的时间t的k阶多项式; ε(t)是观测误差。

在对连续15 d观测数据的分析中,重点是验证数据的完整性、干扰因素以及异常值是否存在。在进行噪声水平分析之前,没有对原始数据进行任何预处理。随机选取蓟县地震台gPhone、DZW和GS-15型重力仪连续15 d的原始重力观测数据,对其使用1.1节中数据处理方法进行处理,以此评估不同型号重力仪在亚地震频段对环境噪声的响应差异。图1分别给出了gPhone、DZW和 GS-15型重力仪最平静的连续15 d的原始重力数据、扣除固体潮和气压改正后的重力数据以及利用截止频率为8 h的高通滤波器滤波后的重力残差功率谱。

图1 3台重力仪连续15 d的原始重力观测值(a)、扣除理论固体潮和大气改正后的重力残差(b)和使用截止频率为8 h的高通滤波器滤波后的重力残差(c)
Fig.1 Original gravity observations from the three gravimeters in a continuous 15-day period(a),gravity residuals after atmospheric corrections and subtracting theoretical solid tide(b),and gravity residualsafter applying a high-pass filter with a cutoff frequency of 8 hours(c)

从图1可以看出,gPhone和GS-15型重力仪的原始重力观测曲线的振幅较为相近,大约为300×10^-8 m/s²,DZW型重力仪记录到的原始重力数据有明显的漂移。扣除固体潮并进行气压改正后,数据质量有所改善,3种型号重力仪的重力残差振幅分别为10×10^-8、250×10^-8 和100×10^-8 m/s²,其中DZW型重力仪重力残差信号振幅最大、漂移最为明显。

本文选取2022年蓟县地震台3种型号重力仪的连续重力观测数据,分别进行在亚地震频段的背景噪声水平计算,发现gPhone、DZW和GS-15型重力仪的固体潮数据平静期分别为7月30日—8月13日、7月4日—18日、7月29日—8月12日。3台仪器平静期较为相近,根据RMS值选取的gPhone、GS-15型重力仪的平静期前后仅相差1 d。3台重力仪最平静15 d的固体潮数据,经过气压改正、高通滤波后得到重力残差,并利用101个点的Parzen窗对其做作平滑处理,得到重力残差的功率谱,结果如图2所示。

图2 蓟县地震台3台重力仪固体潮观测数据平静期重力残差的功率谱密度
Fig.2 Power Spectral Density(PSD)of gravity residuals of three gravimeters at Jixian Station during the calm period

从图2可以看出,在亚地震频段,gPhone型重力仪的噪声水平低于DZW和GS-15型。根据1～6 h频段的平均功率谱,得到gPhone、DZW和GS-15型重力仪的SSNM值分别为3.91、4.48和5.20。

图3给出了25台重力仪在亚地震频段的背景噪声水平,其中超导重力仪结果源自GGP(Global Geodynamics Project)中19个台站(Rosat et al,2004)以及拉萨台和武汉台(张苗苗等,2015)。根据图3可知,25台重力仪中,有16台超导重力仪的背景噪声水平低于新低噪声模型(NLNM)(Peterson,1993),一半以上的超导重力仪的SSNM值小于3。相比拉萨台C032型超导重力仪,武汉台iGrav-007和OSG-065型超导重力仪较新,噪声水平整体偏低。蓟县台3台重力仪的SSNM值均高于NLNM值,其中,gPhone型重力仪噪声水平高于95%的超导重力仪,其SSNM值约为Mo台站的2倍。

将2022年蓟县台3种不同类型重力仪在亚地震频段和地震频段的背景噪声水平(李响等,2024)汇总于表2。从表2可发现,SSNM值和SNM值呈正相关,这与Rosat(2005)得出的SNM值和SSNM值具有线性强相关性的结论一致。因此可以推断,地震频段的背景噪声越大,亚地震频段的背景噪声也越大,台站背景噪声水平就越高,反之亦然。

图3 重力仪在亚地震频段的背景噪声水平
Fig.3 Background noise levels of gravimetersin the seismic sub-frequency band

表2 2022年蓟县台3种类型重力仪在亚地震频段的背景噪声水平
Tab.2 Background noise levels in the seismic sub frequencyband of gravimeters at Jixian Station in 2022

对gPhone、DZW和GS-15型重力仪在潮汐频段的背景噪声水平进行比较分析。首先对2022年的原始观测数据(1 min采样率)使用1.2节中数据处理方法进行预处理,蓟县台3台重力仪预处理后的重力数据和气压数据,如图4所示。

使用Eterna方法对预处理后的数据进行调和分析。采用Tamura(1987)提出的潮位展开表,运用高通滤波剔除仪器漂移的影响,同时对气压变化的影响进行了改正。刘清超等(2016)在南极中山站利用LCR-ET21型弹簧重力仪也进行过类似的研究,表3给出了4种重力仪在不同潮汐频段的噪声水平进行对比分析。

图4 蓟县台重力仪原始数据经过间断、突跳、阶跃修正后的小时采样重力数据和气压数据
Fig.4 Hourly-sampled values of gravity after intermittent correction,abrupt-change correction,and step-change correction at Jixian Station,and barometric recordings

表3 不同重力仪不同潮汐频段的背景噪声水平
Tab.3 Background noise levels of gravimeters in different

从表3可以看出,除去1/4日频段,gPhone型重力仪的噪声水平普遍低于DZW和GS-15型重力仪。DZW型重力仪在所有潮汐频段的噪声水平均低于GS-15型重力仪,特别是在周日频段,差异更为显著,可达5.7倍之多。可能由于南极噪声比较大,gPhone型和DZW型重力仪在所有频段的噪声水平都低于中山站的LCR-ET21型重力仪。

表4列出gPhone、DZW和GS-15型重力仪的观测资料通过调和分析得到的主要潮波(O₁、K₁、M₂、S₂)的潮汐参数,包括振幅因子δ和相位滞后Δφ及其均方根误差。由表4可见,gPhone型重力仪得到的潮汐参数的均方根误差普遍小于其它2种重力仪。由于GS-15型重力仪仪器老化严重,本文着重比较gPhone和DZW型重力仪的潮汐参数。

gPhone和DZW型重力仪的O₁和K₁潮波的差距接近3倍,M₂和S₂潮波的差距达到2倍多。实际上,gPhone型重力仪在调和分析中得到高精度的潮汐参数,如O₁波的精度为0.33‰,K₁波的精度为0.26‰,M₂波的精度为0.24‰,S₂波的精度为0.53‰。

表4 蓟县台3台重力仪调和分析得到的主要潮波潮汐参数
Tab.4 The main tidal parameters of the 3 gravimeters at Jixian Station obtained from the harmonic analysis

表5给出了gPhone、DZW和GS-15型重力仪在地震频段、亚地震频段和潮汐频段的背景噪声水平(Rosat et al,2004; 李响等,2024)。从表5可以看出,在地震频段、亚地震频段和潮汐频段,相较于DZW和GS-15型重力仪,gPhone型重力仪的背景噪声较低。因此,TNM值和SNM值具有一定的相关性,与Rosat等(2004)计算得出19个台站的超导重力仪的SNM值和TNM值具有线性相关性的结论一致。

图5给出了蓟县台3台重力仪和GGP中19个台站超导重力仪(Rosat et al,2004)在潮汐频段背景噪声水平的对比情况。从图5可看出,19台 超导重力仪中,约50%的超导重力仪的背景噪声水平小于5,80%的超导重力仪的背景噪声水平小于5.5。gPhone型重力仪在潮汐频段噪声水平高于部分超导重力仪。

表5 2022年蓟县台3台重力仪在不同频段的背景噪声水平
Tab.5 Background noise levels of the 3 gravimetersin the tidal frequency band in 2022

图5 重力仪在潮汐频段的背景噪声水平
Fig.5 The background noise level of gravimetersin tidal frequency band

选取2022年蓟县地震台3种重力仪的小时采样的重力数据进行平均滤波和降采样处理,以达到分钟的采样频率,并结合该台站的气压观测数据,分别运用Eterna和Venedikov调和方法进行潮汐分析。两种分析均使用了Tamura分波表中的1 200个潮波进行计算,均得到各个潮汐波群的潮汐参数(Tamura,1987)。固体潮参数指标是评估观测数据质量的关键工具,包括Nakai检验百分比、M₂潮汐因子、M₂潮汐因子中误差和相对中误差等。由于蓟县台位于中纬度地区,M₂波在中纬度的振幅最大,所以通过对比观测固体潮M₂波振幅与理论固体潮M₂波振幅的比值及其标准差、M₂潮汐因子与M₂潮汐因子中误差的比值来评价台站数据质量的可靠性。

表6和7给出了分别采用Eterna方法和Venedikov方法对蓟县地震台的重力固体潮观测数据进行调和分析的结果,并提供了详细的波群信息,包括每个波群的起始频率和终止频率、振幅因子及其均方根误差以及相位滞后及其均方根误差。均方根误差是观测值与理论值之间差异的平方的平均值,是衡量观测数据集中趋势的指标,反映了观测误差的平均水平。均方差越小,表示观测数据的波动越小,观测结果越稳定。表6、7显示,所有的波群中理论振幅最大的波群为M₂、O₁、K₁、S₂,而相应的振幅因子中误差远小于其他潮波的中误差,因此潮波振幅越大,中误差越小。

表6 采用Eterna方法对蓟县台重力仪重力数据进行调和分析
Tab.6 Harmonic analysis of gravimeters' recordings at Jixian Station through Eterna method

表7 采用Venedikov方法对蓟县台重力仪数据进行调和分析
Tab.7 Harmonic analysis of gravimeters' recordings at Jixian Station through Venedikov method

Eterna方法采用了最小二乘谐波分析,可精确分离各潮汐波,更适用于长期连续观测。Venedikov方法基于潮汐频带滤波,将数据分成不同频段后进行分析。由于使用的分析方法不同,每种波群的频率范围也略有不同,经过调和分析后,得到的振幅因子的中误差也有差别。本文选取的是一年的数据,对于振幅比较大的M₂、O₁、K₁、S₂,采用Eterna方法得到的振幅因子中误差比较小。