0 引言

地球因太阳、月亮相对位置变化而产生周期性潮汐变化时，固体地球自身产生周期性的地面重力变化，长周波、半日波和日波合并称为重力固体潮(焦新华等，2009)。对重力固体潮的观测和分析，是在获得日月引潮力作用下的地球变化特征后，得到固体地球内部构造特征(刘清超等，2016)，通过计算来提取重力固体潮中包含的地球物理信息以达到获取重力固体潮信号的目的。海潮是海洋表面受到太阳和月亮引力后产生的周期性涨落现象，其与重力固体潮具有同源性，会产生相同的频谱特征，因此在观测的重力固体潮资料中会夹杂海潮信息(张昭栋等，1989)。分析重力固体潮汐残差信号可获得地球内部动力学信息，但其中夹杂的全球和局部海潮负荷对分析结果影响较大，仅采用滤波法很难将其剔除，尤其在沿海地区，海潮负荷对重力潮汐振幅占比更大(杜文成等，2018)。为了更加有效地利用观测到的重力资料精密测定地球潮汐常数以及研究地核运动等问题，对重力观测资料进行海潮负荷信号改正尤为重要(孙和平等，2002)。

Farrell等(1972)首次提出采用海潮潮高和负荷格林函数进行褶积积分计算，求取海潮负荷影响，而后Schwiderski(1980)采用流体动力学差值方法构建较高精度的SCW80海潮模型，被广泛认可并作为标准使用。伴随现代高精度Topex/Poseidon(T/P)卫星测高技术的进一步发展，国内外学者在前人基础上对海潮负荷开展更广泛的研究，通过构建和优化海潮负荷模型，对重力固体潮、地倾斜固体潮的影响进行改正计算，结合地域特点来选取某区域的最优海潮模型，提出了多种高精度海潮模型(Sun，1992；Neumeyer et al，2005；Bos et al，2005；周江存等，2005、2007；周江存，2008；孙和平等，2005a、2005b、2006；Florsch et al，2009；Arnoso et al，2011；Kim et al，2011、2013；Zhou et al，2013；Ducarme et al，2014；杨锦玲等，2016；杜文成等，2018；刘清超等，2016)。本文主要搜集和研究现阶段应用较广泛的FES2004、GOT4.7、EOT11a、DTU10、NAO99b、HAMTIDE11a、TPXO72 atlas和TPXO72等8个海潮模型。

辽宁位于中国东北地区南部，构造上处于郯庐断裂带中北部渤海与辽河平原之间的营潍断裂带上，该地区地质构造复杂，基底介质刚度较高，新生代活动较强，复杂的构造环境使得辽宁地区成为中国大陆东北部最为活跃的地区(翟丽娜等，2019、2020)。大连台gPhone型连续重力观测仪放置于大连市国家基准台观测山洞内，是距海岸线较近的全国重力观测台站，海潮对该台重力固体潮资料的影响较大。为了更加真实地反映海潮影响程度，在对重力数据进行干扰处理后，进行重力固体潮各波群潮汐参数及精度估算，选取适合大连地区的8个海潮模型对主要潮波S₂、M₂、K₂、N₂、K₁、O₁、Q₁和P₁进行海潮负荷改正计算分析。

1 预处理调和分析

自2010年起，大连台gPhone秒采样弹簧重力仪进行观测并连入陆态网络数据库中。本文选取2015年1月1日~2019年8月14日的秒值数据进行分析和研究(图 1、图 2)。进行潮汐参数计算时，首先进行降采样预处理，将秒值降采样为分钟值、分钟值降采样为小时值后，才可进行调和分析。采用国际低潮中心(ICET)降采样滤波器对秒值进行降采样，得到分钟值后，一般采用Tsoft软件(Van Camp et al，2005)对降采样后的数据进行预处理。数据预处理时，需先移除理论合成潮信息，而后扣除大气效应，剔除地震、尖峰、仪器掉格、突跳等错误数据(徐伟民等，2012；韦进，2012)。采用插值方法对仪器故障、断电等数据进行间断处理，将结果加上已扣除大气效应和理论合成潮的信息，最终得到重力潮汐数据。由于数据采集器原件故障，导致降采样处理后的重力数据在2018年9月15日~10月31日期间出现了数据缺失。采用Etema3.30(Wenzel，1996)软件对预处理的重力固体潮整点值做调和分析，计算获得N₂、K₂、K₁、M₂、S₂、O₁、Q₁和P₁波精度较高的参数(振幅因子、相位滞后和误差估算)，见表 1。其中，振幅较大的有半日波S₂和M₂以及周日波K₁和O₁，振幅分别为S₂波28.9938μGal、M₂波48.1454μGal、K₁波42.3419μGal、O₁波29.3428μGal，8个主波振幅因子标准差在3.2%以内，精度较高。

2 海潮负荷计算及分析 2.1 海潮模型

近年来，全球现代化海潮模型有同化模型和经验模型2种，经验模型是对卫星测高数据进行潮汐信息提取而成，同化模型是在流体动力学模型基础上对各类观测数据同化后形成。本文选取8个海潮模型进行分析，虽然分辨率和精度逐步提高，但对于近海地形状况，海潮梯度一般变化较大。对大连连续重力数据进行分析时，选用的海潮模型具体信息如表 2所示，其中，EOT系列海潮模型由德国大地测量研究所提供，EOT11a模型以FES2004模型为基础，选取T/P、ERS-2、Jason1/2和Envist等卫星的测高数据，分辨率为0.125°×0.125°；TPXO系列模型由美国俄勒岗州立大学提出，分辨率为0.25°×0.25°；DTU10模型由丹麦技术大学在法国潮汐小组提出的FES2004模型基础上建立，分辨率为0.125°×0.125°；HAMTIDE系列模型由德国汉堡大学提出，采用最小二乘法在广义反演基础上降低模型与数据的不确定性，分辨率可达0.125°×0.125°；TPXO72 atlas模型是对T/P、ERS-2、Jason-1等3个卫星测高数据做调和分析，结合沿海地区分辨率进行结算，精度为0.25°×0.25°；NAO99b模型是日本天文台应用T/P卫星数据，结合日本验潮站资料构建的区域模型，精度为0.25°×0.25°。

2.2 海潮的负荷效应计算

根据Farrell等(1972)提出的理论，结合Agnew(1997)的积分格林函数计算法，将海潮与重力格林函数做褶积积分计算(刘清超等，2016；杨锦玲等，2016；杜文成等，2018)，具体计算公式为

$ \begin{aligned} L(\theta, \lambda, t) &=a^{2} \iint \rho H\left(\theta^{\prime}, \lambda^{\prime}, t^{\prime}\right) \cdot G(\psi) S(A) \sin \theta^{\prime} \mathrm{d} \theta^{\prime} \mathrm{d} \lambda^{\prime} \\ &=a^{2} \iint \rho H(\psi, A, t) G(\psi) S(A) \sin \psi \mathrm{d} \psi \mathrm{d} A \end{aligned} $

(1)

其中，L(θ，λ，t)表示某时刻某点的海潮负荷量，θ和λ分别为计算点坐标，t为时间；ρ为海水密度；a为地球半径；H(θ′，λ′，t′)为负荷点(θ′，λ′)处的瞬时潮高；G(ψ)表示与角距相关的重力负荷格林函数(该函数与负荷勒夫数相关)。

经过8种海潮模型改正计算后，分别得到O₁、K₁和M₂波受海潮影响的负荷信息，如表 3所示。其中，O₁波最大振幅与最小振幅的差值为0.136μGal、K₁波为0.104μGal、M₂波为0.61μGal，重力海潮负荷的振幅分布频段在2.4~2.8μGal之间，总体上8个海潮模型的差异性较小。

2.3 海潮的负荷改正计算

对重力海潮负荷3个主潮波信息进行矢量叠加计算，以减小海潮负荷影响，并进一步剔除该影响，得到改正后各潮波的潮汐参数。海潮负荷改正后的矢量图如图 3所示，其中矢量A^th表示固体潮理论值，矢量A为固体潮观测值，矢量L为海潮负荷(相位、振幅分别为λ和L)，矢量B表示固体潮观测残差，矢量X表示固体潮剩余残差(相位、振幅分别为x和X)。α、β、λ为重力海潮负荷3个主潮波的矢量角。

潮波频率为σ，B、X及潮汐参数δ^c(Δφ^c)表示为

$ \bar B = B(\sigma){{\rm{e}}^{i\beta (\sigma)}} = {{\bar A}^{{\rm{th}}}}(\sigma)\left[ {\delta (\sigma){{\rm{e}}^{i\Delta \varphi (\sigma)}} - {\delta ^{{\rm{th}}}}(\sigma){{\rm{e}}^{i\Delta {\varphi ^{{\rm{th}}}}(\sigma)}}} \right] $

(2)

$ \bar X(\sigma) = X(\sigma){{\rm{e}}^{i\chi (\sigma)}} = \bar B(\sigma) - \bar L(\sigma) = B(\sigma){{\rm{e}}^{i\beta (\sigma)}} - L(\sigma){{\rm{e}}^{i\lambda (\sigma)}} $

(3)

$ {\delta ^{\rm{c}}}(\sigma){{\rm{e}}^{i\Delta {\varphi ^{\rm{c}}}(\sigma)}} = \delta (\sigma){{\rm{e}}^{i\Delta \varphi (\sigma)}} - L(\sigma){{\rm{e}}^{i\lambda (\sigma)}}/{A^{{\rm{th}}}}(\sigma) $

(4)

2.4 海潮负荷效应计算结果分析

根据海潮负荷理论的计算方法和模型分析，利用8个全球海潮模型对大连台重力数据进行海潮负荷改正计算，发现固体潮剩余残差与观测残差相比减小较多。固体潮观测残差经海潮模型改正后的剩余残差，如表 4所示，计算中主要以O₁、K₁和M₂波为例，分析发现上述3种半日波的观测残差值分别为0.499μGal、0.701μGal和1.825μGal，经海潮模型改正后其剩余残差数值范围在0.153~0.203μGal、0.188~0.388μGal和1.088~1.367μGal之间。残差振幅大于理论值，但总体明显减小，初步估计由未进行区域海潮模型的改正计算以及环渤海部分区域环境因素、大气负荷效应和台站区域地质构造等因素引起，后续工作将重点进行渤海区域海潮模型计算的对比分析。

利用O₁、K₁和M₂波的振幅因子来综合评价负荷改正效果，对主要潮波负荷改的正有效性进行统计，如表 5和图 4所示。负荷改正有效性根据孙和平等(2006)给出的由观测残差矢量振幅减去最终残差振幅、再除以观测残差振幅而获得。对海潮负荷改正效果做综合评价发现，应用8个海潮模型进行海潮改正后O₁波的负荷改正有效性为43%~45%，M₂波为48%~53%，K₁波为43%~48%。在大连台重力海潮模型的改正中，结合表 5及图 4的计算结果分析发现，对K₁波和O₁波海潮改正效果最好的为TPXO72 atlas模型，而M₂波潮波改正效果最好的为TPXO72模型。不同模型剩余残差值仍较大，一方面是由于大连海湾深入渤海腹地位置的特殊性，该地海潮影响与其他地区具有较大的差异性，使其模型误差较大；另一方面，重力仪器运行时间长，仪器零漂较大，考虑其为大气负荷影响。在重力固体潮汐数据预处理时，将大气负荷、验潮站信息和近海海潮模型加入修正会得到更好的计算效果。

3 结论

根据8个全球海潮模型，对大连台重力数据受海潮负荷影响大小进行海潮负荷改正计算，基于计算结果，对比海潮模型在大连台的计算效果，得出如下结论：

(1) 根据8个全球海潮模型计算出O₁、K₁和M₂波的重力海潮负荷振幅分布频段在2.4~2.8μGal之间。对主潮波的海潮负荷振幅数值分析发现，O₁波和K₁波的振幅明显小于M₂波，在经过海潮模型改正后，3个主潮波的负荷振幅及相位均有不同程度的减小。

(2) 对负荷改正效果进行综合评价，3个主潮波的负荷改正有效性为43%~53%，K₁波和O₁波海潮改正效果最好的为TPXO72 atlas模型，而M₂波潮波改正效果最好的为TPXO72模型。总体来说，不同模型剩余残差值仍较大，可能与大连海湾位置的特殊性、近海模型存在误差等有关。如采用近海高精度海潮模型进行修正计算，并在预处理方法中尝试新的改进，会得到更好的计算效果。

致谢: 感谢匿名审稿专家对论文进行了细致的审阅并提出宝贵的审稿意见。