摘要：制约时变卫星重力场模型应用效能的一个重要问题是球谐系数中的噪声影响，主要表现为在地球表层质量场信号恢复中受到南北方向的条带和随机噪声。为了有效抑制南北方向条带噪声，本文提出了一种多域组合GRACE去条带滤波器。该滤波器综合S&W去相关滤波器和SSAS滤波器的优

本文内容来源于《测绘学报》2024年第11期（审图号GS京(2024)2421号）

一种顾及地理纬度空间互补的多域组合GRACE去条带滤波方法

1.中国地质大学（武汉）地理与信息工程学院，湖北 武汉 430074

2.

3.

4.

5.

基金项目

摘要

制约时变卫星重力场模型应用效能的一个重要问题是球谐系数中的噪声影响，主要表现为在地球表层质量场信号恢复中受到南北方向的条带和随机噪声。为了有效抑制南北方向条带噪声，本文提出了一种多域组合GRACE去条带滤波器。该滤波器综合S&W去相关滤波器和SSAS滤波器的优势，采用最小二乘法和多通道奇异谱分析，兼顾频谱域和空间域的信号特征。研究结果表明，相比于传统滤波器，本文提出的组合滤波器能够适应不同纬度的滤波需求，有效消除中高纬度和低纬度的南北方向噪声影响，保留更多地球物理信号，提高观测小尺度信号的可能性。

关键词

第一作者：吴晓辉（2002—），男，硕士生，研究方向为卫星重力数据处理。 E-mail：

本文引用格式

吴晓辉, 吴云龙, 徐国栋, 易爽, 刘素兰, 张豹, 钟玉龙, 王力哲.

WU Xiaohui, WU Yunlong, XU Guodong, YI Shuang, LIU Sulan, ZHANG Bao, ZHONG Yulong, WANG Lizhe.

阅读全文

基于快速、高效和几近全球覆盖的优点，重力卫星可确定高精度、高分辨率的静态及时变全球重力场。利用重力GRACE level-2产品提供的球谐系数可以恢复全球质量场，被广泛用于研究地球表层质量变化，包括冰后回弹、陆地水储量变化[1-2]、全球海洋洋流[3]、地震重力变化[4]和极地冰盖[5]等。其结果提高了人类对地球系统物质迁移的认识，但其信号受到严重的南北条带和随机噪声干扰，需要进行精细的滤波预处理。

目前利用GRACE时变重力场进行各类地学应用的前提是尽可能消除噪声影响，其理想状态是最大程度保留信号和抑制噪声，由此陆续衍生了两类滤波器的设计及使用：第一类滤波器致力于从原始信号中分离出噪声，如S&W（Swenson & Whar）去相关滤波器[6]、奇异谱分析滤波器[7]、经验正交函数滤波器[8-9]等，这类滤波器的设计基于噪声的空间或时间特性，通过一定去相关手段从系数中分离噪声；第二类滤波器则通过调整不同阶次球谐系数权重的方式抑制噪声，如高斯滤波器[10]、在高斯滤波器基础上设计的Han滤波器[11]、Fan滤波器[12]、数据自适应Wiener滤波器[13-14]和均方根（root mean square，RMS）滤波器[15]，以及利用先验信息设计的DDK滤波器等[16-17]。

第一类滤波器设计通常采用单一去相关手段处理全局信号，这导致在特定区域取得良好效果的同时，在其他区域可能产生显著的失真。空间滤波器作为第一类滤波器中特殊的存在，通过分析噪声在空间中表现的规律，并采用相应的去相关手段进行滤波。如被广泛应用的去除南北方向条带（north-south striping，NSS）噪声的S&W去相关滤波器[6,18]，仅能有效处理中高纬度地区的条带噪声，而对低纬度地区的作用效果较差；空间域多通道奇异谱分析（multi-channel singular spectrum analysis in the spatial domain，SSAS）滤波器设计在空间域中去除NSS噪声，但NSS噪声随纬度增加而变稀疏的特点，导致其无法很好处理中高纬度的信号。而第二类滤波器采用调整系数权重的方式对噪声进行压制，不加选择地处理信号和噪声的特点，有效抑制噪声的同时，会带走小尺度信号，降低信号的分辨率。同时由于第二类滤波器会修改低阶系数的权重，而低阶系数中仅包含少量噪声，会损失部分大尺度信号。尽管采用第二类滤波器能够有效去除噪声，且在此基础上设计组合重构滤波器使其信号损失降低，有效提高信噪比[19-20]，但仍然无法解决其降低信号分辨率的问题。

总体而言，这些方法在全局信号滤波均存在不足。如何结合两类滤波器的优势处理区域，在消除噪声的情况下尽可能保留信号特征，最终在全球尺度和特定区域达到最佳去条带效果，是进一步提升卫星重力时变模型地学应用效能的关键预处理环节。

基于此，本文设计了一种顾及地理纬度空间互补的多域组合去条带滤波方法，有效克服了以上两类滤波器的不足，并最大程度消去NSS噪声。通过对不同空间滤波器的去噪特点及作用空间域进行分析，分区采用不同的滤波方法，以更好地处理信号。本文组合使用S&W滤波器和SSAS滤波器的去相关原理，构建了一种顾及地理纬度空间互补的多域组合GRACE去条带滤波方法（de-striping in the spatial and spectral domain，DSS），同时在频域和空间域消除噪声。通过模拟信号和真实信号分别测试，与单一空间滤波方法比较，检验了本文提出的改进算法的效果和可靠性。

本文研究工作基于得克萨斯大学空间研究中心（Center for Space Research at the University of Texas，CSR）提供的RL06和RL06.1 GSM月模型产品，时间跨度为GRACE和GRACE-FO任务的2002年4月至2023年5月。一阶项球谐系数使用TN-13结果替换[21]，C20和C30由TN-14结果替换[22-23]，缺失月份数据采用多通道奇异谱分析插值法进行插值[24]。

将每月的球谐系数转变为全球质量场，以等效水高的形式表示

滤波器设计流程如图1所示，针对不同噪声含量的球谐系数和纬度区域区间，采用多项式去相关和多通道奇异谱分析去相关，分别在频谱域和空间域中去除噪声。首先，对各区间球谐系数的噪声含量进行分析，从而提取原始球谐系数中各区间噪声含量的信息，即分析得到噪声含量高和含量低的系数区间和纬度区间。然后，在频谱域进行去相关滤波，结合各区间系数噪声含量信息和球谐系数的标准差数据对原系数进行分区，针对噪声含量高的系数采用强去相关手段，针对噪声含量低的系数采用弱去相关手段。最后，在空间域中进行去相关滤波，针对噪声含量较高的纬度区间进行去相关滤波，在全球尺度上达到最佳去噪效果。

图1 滤波数据处理流程和各步骤所提取噪声的空间域结果

Fig.1 Flowchart of the filtered data processing and the spatial domain results of the noise extracted by each step

首先对球谐系数中各波段中包含的噪声进行分析，采用P3M10去相关滤波器（S&W滤波器的修改版，采用三次多项式拟合，起始去相关次数order=10）[4]，提取各波段的噪声（原始信号与滤波信号的差值），结果如图2所示。根据球谐系数进行分区，球谐系数的次数项为[10，25）、[25，40）、[40，60]的噪声提取结果如图2（a）、（c）、（e）所示，可以看出随着球谐系数次数增加，噪声所在的纬度降低，强度随纬度降低而降低；球谐系数的阶数项为[10，30）、[30，50）、[50，60]的噪声提取结果如图2（b）、（d）、（f）所示，可以看出低阶系数提取的噪声多为大尺度信号，高阶系数包含更多小尺度信号，且大多数噪声包含在阶数项[50，60]这一区间内。由数据分析得知，低阶系数含有很少噪声，随着阶数增加，噪声增加的趋势明显，且这种趋势并不是线性的。

图2 分波段提取噪声结果及使用P3M10和SSAS方法提取的噪声

Fig.2 Noise extracted in sub-bands and noise extracted using P3M10 and SSAS methods

图2（g）、（h）分别展示了P3M10和SSAS滤波器所提取的噪声。P3M10提取的噪声随纬度降低而减少，说明其能够有效处理低次项系数中的噪声，但研究表明高次项系数中包含的噪声更多，即低纬度区域应当包含更多噪声。SSAS与P3M10在空间域表现结果相反，噪声随纬度增加而减少，其提取的噪声主要位于高次系数。

结合图2展现的噪声在系数中的分布区间，应当重点处理高阶/次系数，即对这部分噪声施加更强的滤波手段。因此本文采用了一种差异化的策略，即对低阶系数采用较弱去相关策略，对高阶系数采用强去相关策略，特别对高阶次的系数采用2种去相关滤波手段以最大程度去除噪声。

综上，本文提出重构组合S&W方法和SSAS方法以在全局和区域上达到最优去噪效果。特别对于噪声较多的低纬度地区，采用最小二乘法和多通道奇异谱分析两种方法进行处理，以最大程度去除噪声，保留更多小尺度信号。

本文采用S&W去相关滤波方法作用于频谱域去噪，控制S&W方法去噪强度的参数包括多项式阶数x和参与去相关的参数个数w，现有研究表明，利用这一方法设计的滤波器根据选择的参数不同，滤波结果也各有特点[4,6,25-26]。文献[6,26]对低阶次球谐系数采用高强度去相关滤波处理，导致有效去除噪声的同时带走了部分中高纬度地区的信号；而文献[4,25]对高阶球谐系数与低阶球谐系数采用同一强度的去相关滤波，造成高阶/次噪声无法有效去除，因此需要较强的后处理滤波器处理。两种策略的共同点是都直接或间接导致了小尺度信号的丢失，区别在于由于文献[6]的策略对近对角线及高次系数采用了更强的滤波手段，相较于文献[4]方法有效减少了低纬度的噪声，使剩余噪声分布在30°S—30°N这一区间内。

考虑到标准差较大的系数表明存在更多不确定性和噪声，利用GRACE level-2产品提供的球谐系数的标准差数据，计算不同时间段球谐系数标准差的平均值，分别代表GRACE运行前10年、GRACE运行后半段时间和GRACE-FO运行至今的误差模式，结合图2中所展示信息对球谐系数进行分区，结果如图3所示。Ⅰ区中球谐系数代表大尺度和噪声较少的信号，对这部分系数采用较弱的去相关手段处理，Ⅱ区中球谐系数代表小尺度和噪声较多的信号，采用强去相关滤波器以在最大程度上消去系数中噪声。

图3 Clm的误差模式（按×1013比例缩放）

Fig.3 The error (scaled by ×1013) pattern of Clm

因此，综合3个运行阶段及考虑后续使用该方法进行各类应用使用的数据精度，选择对Ⅰ区施加强度较弱的滤波以保留大尺度特征，Ⅱ区采用强大的滤波手段压制噪声以增加观测小尺度信号的可能性，因此分别对两个区域采用式（2）—式（6）处理。通过图2和图3展示的结果，确定采用双扭线拟合左侧曲线，并加入直线使区域封闭，如式（2）所示

(2)

式中，θ0为曲线角度；m0为曲线取极大值点的球谐系数次数；l0为当曲线取极大值点处球谐系数阶数；ρ为曲线的极径，此处取ρ=42；m0:l0指从m0到l0所有值。

采用扇形曲线拟合右侧曲线，如式（3）所示，右侧曲线选择（55，0）→（45，45）进行拟合

(3)

式中，l1为当球谐系数次数m=0时的取值，此处取l1=55；r取固定值3.5；β为根据曲线另一端点坐标计算得到的结果，此处β=-10/453.5。

考虑到对低次项系数进行去相关滤波会影响高纬度地区的空间域结果，设置起始次数order=10。并利用式（4）提取原始球谐系数中的噪声

利用式（5）对图3中Ⅰ区的噪声进行提取

利用式（6）对图3中Ⅱ区的噪声进行提取

(6)

式中，w为选择的窗口宽度，用于确定所选择的n的范围。此处取阶数p=2，滑动窗口宽度w=5。

经过频谱域去相关处理后大部分噪声被消去，但空间域结果显示低纬度地区存在大量NSS噪声，因此需要引入第2个滤波器处理剩余条带。低纬度区域采用SSAS方法处理剩余条带，SSAS滤波器通过寻找空间域格网之间的相关性去除NSS噪声，识别条带的效果受到噪声含量较少区域的影响。与针对全球的SSAS方案不同，由于经过第一个滤波器处理后NSS噪声主要分布在30°S—30°N区间，因此仅对这部分的格网值进行处理，并将纬度的采样率提高至1°，在最大程度上识别NSS噪声。通过分析信噪比等指标并结合其在空间域的滤波效果，确定采用30阶以上的球谐系数以放大噪声，避免对大尺度信号造成影响。具体步骤如下。

（1）去除各经度序列的平均值，并针对苏门答腊岛附近区域采用掩膜处理防止地震信号被当作噪声消去。

（2）将处理后的序列构造轨迹矩阵Z，利用奇异值分解Z=UAVT，分解为U、A、V3个矩阵。

（3）寻找U矩阵前20列耦合模态，其中相关系数设置为0.9以获得最佳效果，利用耦合模态重构得到噪声矩阵N。

（4）取Z矩阵和N矩阵得到的差值作为新的Z，重复步骤（2）—步骤（3），并记录每一次得到的1 N矩阵。

（5）将步骤（4）得到的所有N矩阵相加并重构，得到结果为中低纬度地区的噪声。

如图1所示，经过两个去相关滤波处理，能够有效处理空间域中大部分NSS噪声，特别对强噪声的低纬度区域采用空间互补去噪的策略，能够减轻额外的后处理滤波器对真实信号的影响。

为了验证本文提出滤波器的有效性和可靠性，采用全球陆地数据同化系统（global land data assimilation system，GLDAS）NOAH V2.1模型数据模拟对2002年12月信号的影响[27]。将土壤水、积雪、冠层等数据转为截断至60阶的球谐系数，并利用等效水高的形式表示，为了便于观测信号的变化，对原信号放大两倍作为模拟信号，结果如图4（a）所示。DDK滤波器的数值越大代表滤波效果越弱，采用经SSAS+DDK4滤波器处理对应月份的GRACE数据得到的噪声用于模拟存在的NSS噪声和随机噪声，结果如图4（b）所示。分别采用P3M10和DSS组合滤波器后处理，并添加SSAS的结果作为对比，并将其结果与原信号作差，结果如图4（c）—（h）所示。

图4 使用模拟信号处理结果对比

Fig.4 Comparison results of using simulated signal processing

与模拟信号对比，可以看到DSS方法有效消除了大部分条带，在陆地上很好地保留了信号；P3M10很大程度上改变了中高纬度的信号的强度，且在中低纬度受到噪声干扰；而单一的SSAS滤波器难以有效处理高纬度地区的噪声，且在各纬度均无法完全消去NSS噪声。DSS在高纬度地区的去噪基于S&W去相关滤波器原理，因此重点分析对低纬度地区信号的影响。通过与原信号作差可知，DSS在低纬度地区的剩余噪声多为随机噪声，仅存在少数NSS噪声。与单一SSAS滤波器的处理结果对比得知，DSS滤波器对低纬度地区处理效果更好。对比图4（d）、（f）、（h）可知，DSS和原信号的差异最小，考虑到添加的噪声中包含了部分信号会对原信号造成干扰，DSS方法较为完整地恢复了原信号的形态，并保留了观测小尺度信号的可能性。

本文选择3个月的模型产品进行方法测试，代表不同季节、不同噪声强度下的去噪效果，并采用SSAS和P3M10滤波器的结果作为对比，结果如图5所示。同时在2002年12月中等噪声强度信号下，采用高斯滤波、Fan滤波和DDK滤波进行单一或组合滤波处理，对处理后空间域结果进行对比，如图6所示。

图5 P3M10、SSAS、DSS方法在2002年12月、2008年5月、2012年7月的空间域结果对比

Fig.5 Comparison of spatial domain results of P3M10, SSAS, and DSS methods in December 2002, May 2008, July 2012

图6 2002年12月不同组合滤波的空间域结果对比

Fig.6 Comparison of results for different filter combinations in the spatial domain in December 2002

2008年5月的结果代表噪声较少时滤波效果（图5（b）），经过P3M10处理后在低纬度仍然含有明显NSS噪声；而SSAS和DSS有效对全局信号进行处理。但对比SSAS和DSS在海洋的处理结果可知，SSAS仍然存在少量噪声；对比两者在亚马孙河等区域的结果可知，DSS滤波器能有效针对低纬度区域信号进行优化。此外，考虑到SSAS的去相关效果随纬度增加而衰减，在中高纬度地区表现的强信号受到噪声干扰的可能性也随之增加。

2012年7月空间域结果显示存在明显的条带噪声，3种方法的处理结果如图5（c）所示。P3M10滤波器在45°S—45°N区间内的信息遭到噪声污染，当进行后续应用时难以保证数据的准确性；SSAS滤波器在非洲和大西洋地区仍然存在明显的NSS噪声，此外各大流域的信号均受到不同程度污染；与之相比，DSS的处理结果并未出现大范围的NSS噪声干扰，也在最大程度上消去了噪声的影响。这意味着即使在信号质量较差的月份，DSS依然能够在很大程度上消除噪声的影响，当后续进行相关应用时，能够得到更为可信的结果。

2002年12月空间域结果代表了噪声强度适中时的滤波结果，一定程度上代表大多数月份的去噪效果。P3M10滤波器对低纬度噪声处理结果较差，无形中增加了中低纬度信号的强度；SSAS在噪声强度适中时仅能处理较为明显的噪声，比较中低纬度海洋区域的条带，其去条带效果甚至劣于DSS在高强度噪声下滤波效果；而经DSS滤波处理后保留的噪声大多为随机噪声，考虑到空间域滤波器难以对随机噪声进行精准识别，DSS滤波器在最大程度上消除了NSS噪声对原信号的影响。

图6为在2002年12月中等噪声强度下，利用300 km高斯滤波器对3种去相关滤波处理结果进行空间平滑，与300 km Fan滤波处理结果进行对比，同时比较不同强度DDK滤波器的处理结果。在采取300 km平滑半径下，组合滤波器和Fan滤波器都有效去除了剩余NSS噪声，但后者削弱了真实信号的强度，如苏门答腊岛的地震信号。比较3种去相关滤波器采用300 km高斯平滑的结果，P3M10处理结果在海洋上仍然存在较为明显的NSS噪声，这一现象很大程度上降低了结果的可靠性，SSAS和DSS处理结果并未出现明显的NSS噪声。图6（f）为SSAS和DSS结果的差值，调整色带范围为[-20，20]，差值信号主要表现为NSS噪声，体现了DSS滤波器特别对低纬度地区信号进行优化的效果。DDK滤波数值越大滤波效果越差，采用DDK4滤波对DSS滤波结果后处理，与单一DDK3、DDK4滤波进行对比，结果之间差值如图6（h）、（j）所示。相比于采用单一DDK4滤波处理，采用DSS+DDK4的结果对NSS噪声有更强的压制作用，与DDK3滤波处理结果之间差值说明DSS+DDK4对低纬度海洋地区噪声有更好的压制作用，同时对陆地的信号有更好的保留。

比较DSS、P3M10和SSAS对全局信号的处理可知，DSS滤波器最大程度上降低NSS噪声对原始信号的干扰，这意味着仅需采用1个较弱的后处理滤波器，即可观测到可用的信息。此外，通过与另外两种方法在低纬度的结果对比可知，P3M10得到的结果受噪声干扰最大，SSAS次之，DSS在中低纬度区域的去噪效果提升明显。相较于使用单一滤波手段，在采用后处理滤波处理剩余噪声时，能够有效避免噪声对真实信号的干扰，为相关地区的研究提供了一种更为有效的手段。

频谱域结果的比较提供了评估不同方法差异的重要信息。图7展示了2002年12月的频谱域结果，代表噪声强度适中时滤波处理后结果。如图7（a）所示，原始系数在高阶和近对角线位置呈现异常大值，这是噪声产生的主要原因，因此考察一个滤波器的效果可以从其对高阶系数的处理结果中看出。

图7 2002年12月的每个系数的振幅和各阶振幅比较

注：C为球谐系数。

Fig.7 Comparison of amplitude per degree and amplitude per coefficient for the December 2002

为了评价3种滤波器的去噪效果，本文采用阶振幅这一概念，各阶系数的振幅大小计算公式如下

(7)

式中，σl为阶振幅。计算结果如图7（b）所示，原始数据集的阶振幅在30阶左右趋势逆转，3种方法都能够在一定程度上降低中高阶强度。P3M10方法无法处理高次项的系数，导致其在高阶处振幅显示出异常大值；SSAS和DSS方法均能压制高次项的噪声，但是SSAS滤波器在高阶处出现小幅度回升，说明其有效处理低纬度条带的同时对中高纬度的去噪效果不足；而DSS滤波器有效地结合最小二乘法和多通道奇异谱分析方法去噪的优势，兼顾对中高纬度和低纬度的去噪，特别对低纬度信号进行组合滤波，有效降低高阶振幅回升的问题。

球谐系数中的噪声随阶数和次数的增加而增多，因此优秀的滤波器应当对高阶系数有更强的去噪效果。由图7（c）、（d）可知，P3M10无法处理近对角线系数的噪声及高次项，导致了低纬度地区的噪声无法消除。而SSAS方法虽然很好地解决了以上问题，但由图7（f）可知其对高阶系数的处理力度不足，特别对接近最大阶数附近的系数，去噪效果有明显的减弱。这也是图7（b）中SSAS在40阶以后阶振幅变化不够明显且在高阶处有小幅度回升的主要原因之一，其主要作用对象是NSS噪声，但并不能完全消除NSS噪声。相比之下，DSS滤波器对高阶系数的去噪效果可由图7（h）看出，其对近对角线系数及高阶系数的去噪效果均有明显提升，凸显了DSS优秀的去噪特性。

由图7（d）、（f）、（h）也可看出3种滤波器的特点。P3M10去噪强度随阶数增加而增强，但对高次系数和近对角线系数并未表现出更高的去噪强度，这与滤波器的特点相关，对边缘系数无法做到有效处理。SSAS有效对近对角线系数的噪声进行处理，但观察图7（f）可知，其对高阶次系数的作用强度较差，即对噪声含量最多的系数效果较差。相比之下，由DSS滤波器结果可以看出随阶数增加，滤波强度随之增强，且对近对角线系数进行有效处理；滤波强度随次数增加而增强，表明对低纬度地区噪声有了更好处理，这与NSS噪声的分布规律相符。尽管对低纬度地区采用高强度的去相关滤波可能会带走部分真实信息，但与采用高强度的后处理进行处理相比，前期对NSS噪声处理有助于保留更多小尺度信号。

为了展示DSS在信号保留及空间分辨率方面的改进，比较P3M10、SSAS、DSS方法处理后的结果，并与CSR Mascon[28]、戈达德太空飞行中心（Goddard Space Flight Center，GSFC）Mascon[29]、喷气推进实验室（Jet Propulsion Laboratory，JPL）Mascon[30]3种Mascon产品进行对比。

为了凸显DSS在中高纬度及低纬度区域的改进，选择3个地区（维多利亚湖、纳尔马达河、新疆某地区），采用DDK6滤波后处理，在不干扰有效信号的情况下同时消除剩余NSS噪声。添加单独采用DDK6处理结果展现DDK6滤波的影响，与Mascon数据对比结果如图8所示。在低纬度的维多利亚湖及其附近区域是受到NSS噪声影响最大的区域，P3M10和SSAS在低纬度区间信号受到NSS噪声影响，在图8中表现出南北方向的质量增加信号。除此之外，P3M10结果中存在部分随机噪声，这一现象仅在P3M10和DDK6的结果中出现，说明未经完全处理的NSS噪声会导致信号增强，从而影响结果准确性。DSS结果中并未出现明显NSS噪声和随机噪声现象，比较3种滤波器处理后结果，说明对低纬度地区采用组合滤波的策略能够对NSS噪声进行有效处理，从而提高所观测小尺度信号的可信度。纳尔马达河及其附近区域位于经过S&W滤波器处理后NSS噪声和有效处理噪声区域交界处。P3M10和SSAS在纳尔马达河附近区域的噪声没有被很好地消除，直接导致标记地区的信号强度增加，海洋上的NSS噪声也很好地验证了这一点。DSS滤波结果中不存在明显的NSS噪声，采用DSS+DDK6这一组合即可有效处理NSS噪声，这是P3M10和SSAS都无法做到的。尽管DSS滤波所得信号的强度明显弱于其他方法的结果，但接近于Mascon结果，且无NSS噪声干扰，因此所得结果的可信度明显高于其他方法。新疆某地区的结果很好地展现了DSS对中高纬度信号的优化。3个Mascon数据均表现出强烈的质量增加信号，但P3M10方法在处理时削弱了这部分信号，而SSAS方法在中高纬度地区的去噪强度较弱，导致该地区信号受到随机噪声的干扰而发生偏移。所得结果DSS对新疆某地区信号的强度和位置结果都与3种Mascon数据表现出相似性，进一步增加了DSS处理结果的可信度。同时注意到除DSS外，其他方法在新疆某地区的结果均未表现出强质量增加信号，这一现象表明对NSS噪声进行有效去除的重要性，对NSS信号进行有效去除有助于恢复信号的位置和真实强度。

图8 不同滤波方法和Mascon数据在维多利亚湖、纳尔马达河和新疆某地区的空间域结果对比

Fig.8 Comparison of spatial domain results of different filtering methods and Mascon data in Victoria Nyanza, Narmada river and a region in Xinjiang

3种滤波方法在指定点质量变化时间序列如图9所示，选用Mascon数据和DDK6滤波结果进行对比。CSR Mascon和JPL Mascon数据致力于恢复信号的强度，而对信号的具体位置有所偏移，如图8中纳尔马达河地区和新疆某地区的结果。SSAS方法在高纬度地区受到噪声干扰比较明显（新疆某地区），而P3M10方法在中低纬度地区表现出更高的信号强度（维多利亚湖、纳尔马达河），这都与直接或间接受到NSS噪声的影响有关。DSS方法在3个地区均表现出其优良特性，特别是对于一些小尺度信号的保留，表现更好。

图9 不同滤波方法和Mascon数据在指定点质量变化的时间序列

Fig.9 Time series of mass changes at specified points for different filtering methods and Mascon data

由于GSFC Mascon恢复信号的位置与利用球谐系数恢复信号的位置相近，且3种方法的时间序列与GSFC Mascon数据结果相近，因此在图9对应图例处标注3种方法与GSFC Mascon数据差值的均方根。本文分别计算了与GSFC Mascon数据差值的均方根，以及与3种Mascon产品数据均值差值的均方根，结果见表1和表2。由表1可知，DSS在维多利亚湖计算的均方根值为6.076 3，相较于另外两种方法有显著提升。在纳尔马达河和新疆某地区分别为3.397 3和4.302 3，特别是不采用DDK6后处理结果的均方根值仍然表现优秀，这一结果说明DSS对低纬度地区的滤波效果优于传统的单一滤波策略。由表2也得到相似的结果，注意到DSS在维多利亚湖的结果并不是最低的，因此添加3种Mascon方法与其均值的均方根结果，从中可以看出3种Mascon结果在低纬度地区的信号差异较大。由图9可知，CSR Mascon与JPL Mascon结果有较为强烈的波动，因此对于该地区，本文选用GSFC Mascon结果为标准，DSS更为接近GSFC Mascon结果。这一结果说明对低纬度区域采用组合滤波的方案能够提升所得数据的准确性，受到NSS噪声影响的可能性降低。

表1各滤波方案结果与GSFC Mascon结果差值的均方根

Tab.1

新窗口打开| 下载CSV

表2各滤波方案结果与3种Mascon均值结果差值的均方根（添加Mascon结果作为对照）

Tab.2

新窗口打开| 下载CSV

综上，DSS滤波器相较于P3M10对各纬度信号强度进行了优化调整，增加在中高纬度信号的强度，降低在低纬度信号的强度，使其更接近Mascon结果。相比于单一SSAS滤波器结果，DSS特别优化了对中低纬度噪声的滤波效果，有效提高了观测小尺度信号结果的可靠性；在中高纬度区域有效去除NSS噪声，减少了其对真实信号的干扰。

为了评价不同滤波方法对信号保留的优劣性，利用文献[15]的RMS_Ratio公式计算各方法处理后的信噪比SNR，并采用10×lg对原式进行缩放以识别方法之间的差异，如式（8）所示

(8)

一般情况下为防止陆地信号泄露，设置300 km缓冲区，MassLand和MassOcean分别指陆地和海洋的质量，ErrorLand和ErrorOcean为陆地和海洋的误差，一般认为陆地上误差水平和海洋上相似，所以ErrorLand≈ErrorOcean。

对3种方法采用DDK6滤波后处理以消除剩余NSS噪声，并利用式（8）计算在每一个纬度上的信噪比随时间变化的结果（图10）。通过每个纬度的信噪比图像对比可以得出各部分纬度地区信号对全局信噪比影响的信息，如图10（d）、（f）所示，注意到并非每个纬度区间都存在陆地和海洋，因此在图中存在空白值，如60°S及90°N附近地区。可以看出P3M10在低纬度地区的信噪比相比于其他两种方法有明显的降低，在采用DDK6去除剩余NSS噪声情况下在低纬度地区的信噪比仍然小于DSS；SSAS方法在高纬度地区信噪比较高，有较大原因是其并未处理陆地的噪声，而在高纬度地区陆地面积较大的区域RMS（MassLand）的结果偏大，而在中纬度及低纬度地区的信噪比均低于DSS。与前两种方法相比，DSS滤波器结合了二者的优点，通过比较图10（a）、（c）、（e）可知，DSS有效提高了中纬度地区和中低纬度地区的信噪比。

图10 3种方法在各纬度信噪比的时间序列及经DDK6处理后3种方法全局信噪比随时间变化的时间序列

Fig.10 Time series of the signal-to-noise ratios (SNR) of the three methods at each latitude, and the time series of the global SNR of the three methods after processing by DDK6 over time

图10（b）为3种方法的全局信噪比随时间变化结果，图例中标注数值为3种方法的全局信噪比随时间变化的平均值。可以看出在绝大多数月份，DSS方法的信噪比是3种方法中最高的，平均信噪比达到6.241 7。与P3M10结果相比，在2010—2018年这个时间段，DSS的全局信噪比有较为明显的提升，可归因于对低纬度地区的信号质量的优化。与SSAS结果相比，DSS在2019—2023年全局信噪比有所提升，通过各纬度信噪比结果可知DSS优化了中低纬度地区信噪比。因此，相比另外两种滤波手段，DSS滤波器在很大程度上优化了除高纬度地区之外的信噪比。

噪声的去除作为GRACE球谐产品处理中一个不可避免的问题，选用不同的去噪方法对小尺度信号的保留程度有重大影响。本文提出一种作用于空间域和频谱域的空间互补的多域组合GRACE去条带滤波器DSS，滤波器设计基于最小二乘去相关法和多通道奇异谱分析原理，分别作用于频谱域和空间域系数，有效去除其中存在的NSS噪声。通过结合信噪比等指标，并分析选择不同滤波参数下的空间域结果以选择最优参数，通过多种试验结果验证了DSS方法的有效性，并得出以下结论。

（1）本文提出的空间互补组合滤波法，相较于传统的单一空间滤波法，更有效处理全局NSS噪声。模拟信号和真实信号下的空间域结果比较展示了不同滤波方法的全局去噪效果，区域结果的比较展现了不同滤波器对小尺度信号的保留效果，结果表明组合滤波器无论在全球尺度上还是对小尺度信号的保留信号效果都表现良好。

（2）对比3种滤波器的信噪比和频谱域结果可知，DSS滤波器在大多数月份拥有最大的信噪比，分析不同纬度的信噪比可知中低纬度地区的信噪比有较为明显的提升。每个系数振幅结果表明DSS滤波器去噪强度随球谐系数阶数增加而增强，阶振幅结果表明DSS滤波器结合S&W滤波器和SSAS滤波器特点，有效处理高阶系数中包含的噪声，从而保留更多有效信息。

DSS滤波器采用空间互补的思路，通过利用不同滤波器的特性，结合同类型滤波器的优势区间，有效对全局信号进行处理。但鉴于S&W滤波器会在高纬度区域（如格陵兰岛）产生一定失真现象，因此DSS滤波器主要用于中纬度及中低纬度的信号处理。在低纬度地区，相较于采用单一滤波器，DSS滤波器的去噪效果有明显改进，为获得该纬度地区更准确的数据提供了一种可用的手段，有望在各个领域的后续应用发挥重要作用。

来源：测绘学报