重力恢复与气候实验（gravity recovery and climate experiment，GRACE）卫星时变重力数据虽然已经被广泛应用于南极冰盖质量变化的估算^[1-4]，但不同估算结果间的差异较大，导致差异的主要原因为冰川均衡调整（glacial isostatic adjustment，GIA）模型不同和泄漏误差改正方法不同。由于人类目前对冰盖历史变化和地幔粘滞度结构认识不足，已有的GIA模型存在较大的误差和不确定性^{[5, 6]}。因此，有效改正南极泄漏误差将会极大提高冰盖质量变化估算结果的可靠性。受GRACE卫星轨道误差和高频背景模型误差等因素的影响，时变重力场模型在发布之前要对球谐系数进行截断。此外，后期计算时，仍要对球谐系数进行滤波去噪处理，而球谐系数的截断和滤波均会造成地球物理信号的泄漏^[7]。Chen等^[2]采用全球正向建模恢复法改正全球区域的泄漏误差，该方法已经被广泛应用于南极冰盖质量变化的估算中^{[4, 7]}。全球正向建模恢复法改正南极泄漏误差的效果与初始信号的准确性密切相关，球谐系数的截断和滤波处理导致南极周边海洋信号内泄漏至相邻的南极陆地，造成相应区域的初始信号不准确，前期研究忽略了该因素对南极陆地泄漏误差改正的影响，影响了冰盖质量变化的估算。

本文利用CSR Mascon数据进行模拟实验，讨论周边海洋信号内泄漏对南极泄漏误差改正的影响；根据实验结果确定采用“流域”函数分离海洋信号，在此基础上，利用全球正向建模恢复法对CSR RL05 GSM（glacial systems model）数据进行处理，并将处理结果与CSR发布的Mascon数据计算结果进行对比。

1 数据与方法

1）GRACE数据。本文采用的数据为CSR发布的2003-01—2013-12总计124个月的RL05 GSM数据（最大阶数为60阶）和Mascon格网数据，其中，2003年6月、2011年1月、2011年6月、2012年5月、2012年10月、2013年3月、2013年8月、2013年9月数据缺失。在对GSM数据进行预处理时，加入了一阶项系数^[8]，用通过卫星激光测距（satellite laser ranging，SLR）解算的C₂₀代替GSM系数中的C₂₀项^[9]，并采用GW13模型扣除冰后回弹影响^[10]。

2）计算方法。地表质量密度变化可由时变重力场模型球谐系数计算得到^[11]，计算公式如下：

$ \Delta \sigma(\theta, \lambda)=\frac{a \rho_{a}}{3} \sum\limits_{l=0}^{\infty} \sum\limits_{m=0}^{l} \bar{P}_{l m}(\cos \theta) \frac{2 l+1}{1+k_{l}} \times\left(\Delta C_{l m} \cos (m \lambda)+\Delta S_{l m} \sin (m \lambda)\right) $

(1)

式中，a为地球半径；ρ_a = 5 517 kg·m^-3，为地球平均密度；l和m分别为球谐系数的阶和次；P_lm(cosθ)为完全规格化缔合勒让德函数；ΔC_lm和ΔS_lm为时变重力场模型球谐系数的变化值；k_l为阶数l对应的负荷勒夫数；θ和λ分别为地心余纬和地心经度。地表质量密度变化通常以等效水高的形式表示：

$ \Delta h(\theta, \lambda)=\Delta \sigma(\theta, \lambda) / \rho_{w} $

(2)

式中，ρ_w = 1 000 kg·m^-3，表示水的密度。

若地表质量密度变化已知，则可反推其对应球谐系数的变化值^[12]：

$ \left\{\begin{array}{c} \Delta C_{l m} \\ \Delta S_{l m} \end{array}\right\}=\frac{a^{2}}{M_{e}} \frac{k_{l}+1}{2 l+1} \sum\limits_{i=1}^{360} \sum\limits_{j=1}^{720} \Delta \sigma\left(\theta_{i}, \lambda_{j}\right) \bar{P}_{l m}\left(\cos \theta_{i}\right)\left\{\begin{array}{l} \cos \left(m \lambda_{i}\right) \\ \sin \left(m \lambda_{i}\right) \end{array}\right\} \sin \left(\theta_{i}\right) \Delta \theta \Delta \lambda $

(3)

式中，M_e = 5.97.19 × 10²⁴ kg；Δθ = π/360；Δλ = π/360（0.5° × 0.5°格网）。

3）南极区域划分。本文按照如下方法将南极区域进行划分：从89.5°S开始，以1°为间隔，用i（i = 1，2，3，···，n）向北计算圆圈数，第i个圆圈上均匀分布着6 × i个格网点^[13]。结合ICESat流域数据，将整个南极分为南极半岛（图 1中的红色区域）、西南极（图 1中的蓝色区域）和东南极（图 1中的绿色区域）。

2 数值模拟实验 2.1 实验方法

全球正向建模恢复法步骤如下：

1）由时变重力场模型滤波信号M_O合成迭代模型M_T。保持滤波信号中陆地数据不变，海洋信号等于负的陆地信号加权值之和的平均值。

2）将迭代模型球谐展开并截断至60阶，采用与步骤1）中相同的滤波方法得到迭代模型滤波信号M_P，计算M_O与M_P的差值ΔM。

3）将ΔM乘上加速因子k(k = 1.2)后，加入迭代模型M_T中：M_T = M_T + ΔM × k。

4）更新的迭代模型即为步骤1）的输入值，重复步骤1）~步骤3），当ΔM小于迭代设计的阈值时，迭代结束，此时M T即为滤波信号经泄漏误差改正后的结果。

根据上述实验流程，可得到经泄漏误差改正后的月平均质量变化的时间序列，对其进行最小二乘拟合，可得到南极冰盖质量变化的长期趋势。

2.2 海洋信号内泄漏对泄漏误差改正的影响

本文以2003—2013年的CSR Mascon数据为模拟数据，按如下方法讨论海洋信号内泄漏对南极泄漏误差改正的影响：

1）以全球Mascon数据为模拟数据，将模拟数据球谐展开，并截断至60阶，经500 km的高斯滤波后转为等效水高数据，按§2.1中的实验流程恢复泄漏信号。

2）Mascon数据在南极区域至少存在120 km的信号泄漏^[14]，为充分考虑自身泄漏对恢复结果的影响，本文以扣除海岸线外200 km内的海洋信号后的Mascon数据为模拟数据，采用与步骤1）中相同的处理，并将结果与步骤1）中的恢复结果进行对比。

图 2为滤波信号、恢复信号和模拟信号长期变化趋势的空间分布图，相应的质量变化速率统计信息见表 1。在空间分布上，扣除海洋信号后的恢复信号与模拟信号一致性较好；海洋信号对东南极影响显著，造成Coats Land、Queen Maud Land、Enderby Land、Kemp Land等地区恢复信号失真。在质量变化量级上，海洋信号影响下的恢复信号质量变化速率为-171.36 Gt/a，与模拟信号相差9.81 Gt/a，扣除海洋信号后的恢复信号与模拟信号相差2.33 Gt/a，且各分区与模拟信号的差异均更小。

造成上述差异的主要原因如下：东南极区域外海洋质量变化信号明显，截断和滤波造成该区域海洋信号泄漏至相邻的南极陆地，导致初始信号在东南极不准确。因此，该区域的恢复信号与模拟信号差异明显。

综上所述，全球正向建模恢复法虽能有效减小南极泄漏误差，但南极周边海洋信号内泄漏会造成东南极的Coats Land、Queen Maud Land、Enderby Land及Kemp Land等地区恢复信号失真，其对全球正向建模恢复法改正南极泄漏误差的影响不可忽略。

3 CSR RL05 GSM数据反演南极冰盖质量变化

本文选取2003—2013年的CSR RL05 GSM数据，采用滤波半径r=500 km的高斯滤波抑制模型条带误差及高频噪声，利用全球正向建模恢复法改正南极泄漏误差。图 3展示了Mascon、GSM滤波信号和全球正向建模恢复信号的长期变化趋势，相应的质量变化速率统计信息如表 2所示。在空间分布上，恢复信号在南极半岛和西南极区域与Mascon结果较为一致，而受东南极区域外海洋信号内泄漏影响较大，恢复信号在该区域与Mascon结果差异明显，且主要体现在Coats Land、Queen Maud Land、Enderby Land和Kemp Land等地区。在质量变化量级上，恢复信号与Mascon结果在南极整体上相差8.2 Gt/a。

上述结论与模拟实验吻合，因此，为准确估计南极冰盖质量变化，在迭代恢复前应扣除海洋信号内泄漏的影响。但现阶段南极周边海洋实测数据缺乏，无法利用外部数据消除其影响，为此本文引入“流域”函数分离海洋信号，以减小海洋信号内泄漏对泄漏误差改正的影响，其推导公式如下^[13]：

$ \Delta h_{\text {ocean }}(\theta, \lambda)=C(\theta, \lambda) \Delta h(\theta, \lambda) $

(4)

当计算格网点位于海洋区域时，C(θ，λ)= 1；位于陆地区域时，C(θ，λ)= 0。对海洋质量变化数据进行球谐展开得到其对应的球谐系数，从GSM系数中扣除，利用剩余球谐系数计算南极冰盖质量变化。要注意的是，GSM数据在南极边界线附近也存在一定范围的信号泄漏，为避免南极泄漏信号被当做海洋信号而被扣除，本文将南极及其海岸线外100 km范围内的信号视为南极陆地信号。

图 4为扣除海洋信号后的全球正向建模恢复结果，相应的质量变化速率信息如表 3所示。由图 4和表 3可知，扣除海洋信号后的恢复信号与Mascon数据在空间分布上基本一致；全球正向建模恢复信号的长期变化趋势为-180.66 Gt/a，与Mascon数据的估算结果仅相差0.51 Gt/a，且各分区与Mascon数据差异均小于1.50 Gt/a。实验结果表明，在无法利用外部数据扣除海洋内泄漏信号时，采用“流域”函数分离海洋信号可有效减小海洋信号内泄漏对全球正向建模恢复法改正南极泄漏误差的影响。

4 结束语

全球正向建模恢复法已被广泛应用于南极泄漏误差改正，但前期研究忽略了周边海洋信号内泄漏对改正结果的影响。本文利用CSR Mascon数据进行数值模拟实验。实验结果表明，海洋信号内泄漏造成东南极Coats Land、Queen Maud Land、Enderby Land和Kemp Land等沿海地区反演结果不准确。本文引入“流域”函数分离海洋信号，并在此基础上利用全球正向建模恢复法处理CSR RL05 GSM数据，估算南极冰盖质量变化，其质量变化的长期趋势在南极整体及各分区上与Mascon数据差异均较小，且空间分布一致性较好。模拟实验和真实数据反演结果均表明，在利用全球正向建模恢复法改正南极泄漏误差前，应排除周边海洋信号内泄漏的影响，才能确保反演结果在质量变化量级和空间分布上的准确性。