時變重力場反映地球系統(tǒng)內(nèi)物質(zhì)的分布、運動和變化，是研究地球系統(tǒng)內(nèi)物質(zhì)運動狀態(tài)及其動力學(xué)機(jī)制的重要約束信息，可為解決人類面臨的資源、環(huán)境和災(zāi)害等緊迫問題提供重要的基礎(chǔ)地球物理信息。近十幾年來，GRACE(Gravity Recovery And Climate Experiment)衛(wèi)星已成功獲取全球均勻覆蓋的時變重力觀測數(shù)據(jù)，使地球時變重力場模型在分辨率和精度上實現(xiàn)了重大提升。GRACE月時變重力場在大空間尺度上定量揭示了全球環(huán)境變化，例如：陸地水儲量變化、冰川消融、海洋變化、冰后回彈以及大地震等導(dǎo)致的地球表層質(zhì)量重新分布。

GRACE重力衛(wèi)星具有可重復(fù)、大范圍的觀測優(yōu)勢，但僅對重力場信號的中長波長部分較為敏感，對短波長部分敏感性較弱，不利于恢復(fù)較高階的時變地球重力場模型。GRACE的實際空間分辨率約400 km，其數(shù)據(jù)存在空間域的南北向條帶等噪聲需做濾波處理，會進(jìn)一步降低空間分辨率，難以滿足環(huán)境保護(hù)、資源監(jiān)測以及災(zāi)害預(yù)測等研究對高精度高分辨率時變重力場的需求。

地面重力時變觀測在地表進(jìn)行，通過定期同點位復(fù)測來獲得重力場隨時間的變化。與GRACE重力衛(wèi)星相比，地面重力測量的優(yōu)點在于重力信號衰減幅度小，觀測結(jié)果的短波長部分更加精確，空間分辨率也更高，適于分析區(qū)域性、近地表的物質(zhì)運動和遷移現(xiàn)象。但地面重力觀測易受地理和環(huán)境等因素制約：一方面，很多區(qū)域（如南極、沙漠和冰川等）難以施測，存在數(shù)據(jù)空白區(qū)，無法像GRACE衛(wèi)星觀測那樣實現(xiàn)數(shù)據(jù)全球覆蓋；另一方面，在可施測區(qū)域內(nèi)，測點空間位置也難以實現(xiàn)均勻分布。直接利用地面重力觀測數(shù)據(jù)確定區(qū)域重力場模型時，由于數(shù)據(jù)僅區(qū)域覆蓋且空間不均勻分布，無法滿足球諧分析這樣的經(jīng)典重力場恢復(fù)理論對數(shù)據(jù)的要求（基函數(shù)在區(qū)域內(nèi)部不滿足正交性）。因此，需借助局部重力場恢復(fù)理論和方法。

目前應(yīng)用較為廣泛的局部重力場恢復(fù)方法包括最小二乘配置法、球冠諧分析法、Mascon法和特殊基函數(shù)法等。Slepian局部譜分析方法即為一種特殊的基函數(shù)法，最早由D. Slepian提出，旨在解決一維連續(xù)情況下的時域和頻域能量集中問題，隨后被不斷發(fā)展并引入到重力場研究領(lǐng)域中(Simons and Dahlen, 2006; Simons et al., 2006; Simons, 2010)。Slepian方法通過構(gòu)建在研究區(qū)域內(nèi)部的正交基函數(shù)，可使信號能量集中在研究區(qū)域內(nèi)部，適用于表達(dá)區(qū)域尺度的重力場變化，并且Slepian展開與相同階次的經(jīng)典球諧展開具有同樣的意義(Simons, 2010; 陳石等, 2017)。

本文對利用Slepian局部譜分析恢復(fù)局部重力場的方法進(jìn)行了研究，并將Slepian方法應(yīng)用到我國華北地區(qū)，確定了華北2011–2013年尺度的120階局部重力場時變模型，并與研究區(qū)域內(nèi)GRACE衛(wèi)星重力估計結(jié)果進(jìn)行對比分析。

Slepian方法的基本原理及恢復(fù)局部重力場模型的流程詳見韓建成等(2021)。Slepian方法是一種具有稀疏(sparse)特性的時空局部譜方法。圖1直觀反映了華北地區(qū)Slepian函數(shù)分別展開到90、120和180階時聚集因子的稀疏分布情況。對于120階Slepian展開，聚集因子總數(shù)量為14641個，其中大于0.001的共41個，僅占總數(shù)量的0.3%，剩余均接近0，分布具有明顯的稀疏性。

圖1 華北地區(qū)對應(yīng)90、120和180階Slepian展開的聚集因子分布

圖2進(jìn)一步給出了華北地區(qū)Slepian展開到90、120和180階時基函數(shù)的空間分布情況，其中洋紅色曲線標(biāo)記重力數(shù)據(jù)覆蓋范圍的邊界。由圖2可以看出，對于每組展開，當(dāng)聚集因子較接近1時，所對應(yīng)基函數(shù)基本都集中在研究區(qū)域內(nèi)部，能量很少泄漏到區(qū)域外部；當(dāng)聚集因子很小時（如圖2第3行中接近0.001的情況），所對應(yīng)基函數(shù)基本集中在研究區(qū)域外部，在研究區(qū)域內(nèi)部幾乎沒有貢獻(xiàn)。

圖2 華北地區(qū)3組Slepian展開的基函數(shù)空間分布，分別對應(yīng)90階（左邊2列）、120階（中間2列）和180階（右邊2列）展開. 對于每組展開，第1和第2行為聚集因子最大的前4個基函數(shù)，第3行為聚集因子接近0.001的基函數(shù). 洋紅色線為研究區(qū)域邊界，黑色線為陸海邊界

本文采用了華北地區(qū)重力測網(wǎng)2011–2013年間的實際重力觀測成果，一共6期觀測數(shù)據(jù)（每年上下半年各測1期）。重力網(wǎng)平差采用基于貝葉斯平差模型的pyBACGS軟件完成(Chen et al., 2019)，該方法可以較好地解決多臺相對重力儀器聯(lián)合平差的權(quán)系數(shù)優(yōu)化、非線性漂移估計和儀器格值的最優(yōu)化估計等問題，6期觀測平差結(jié)果的點值精度均優(yōu)于10×10-8 m/s2. 圖3給出了6期數(shù)據(jù)覆蓋的重復(fù)測點位置的空間分布，共有364個。

圖3 本文研究區(qū)域及重復(fù)測點位置分布. 洋紅色線為研究區(qū)域邊界，綠色十字代表重復(fù)測點位置，灰色線為省界，黑色線為陸海邊界

在研究區(qū)域內(nèi)部364個重復(fù)測點處，將相鄰兩年對應(yīng)測期內(nèi)的重力觀測做差（后一年觀測減去前一年）獲得1年尺度重力變化。所有比較都限定在相對應(yīng)的測期內(nèi)，以減弱水文變化等季節(jié)性信號的影響。最終得到4組1年尺度重力變化，即基于第1期觀測的2013.04–2012.03、2012.03–2011.03（“–”號表示前后兩組數(shù)據(jù)做差），以及基于第2期觀測的2013.08–2012.08和2012.08–2011.08. 獲得研究區(qū)域重復(fù)測點處的重力變化后，接下來應(yīng)用Slepian方法恢復(fù)華北地區(qū)120階時變重力場模型（最佳截斷數(shù)取26），空間變化結(jié)果（格網(wǎng)間距為0.1°×0.1°）如圖4所示。最大展開階數(shù)以及Slepian最佳截斷數(shù)的確定，詳見韓建成等(2021)。

除1年尺度時變重力場模型外，本文還解算了2組華北地區(qū)120階的2年尺度重力變化模型，第1組利用基于上半年觀測(2013.04–2011.03)的重力變化數(shù)據(jù)解算，第2組利用基于下半年觀測(2013.08–2011.08)的重力變化數(shù)據(jù)解算，解算結(jié)果如圖5所示。

圖 4 華北地區(qū)由Slepian方法確定的120階1年尺度重力場變化 (2011–2013) 空間分布：(a)-(c) 基于相鄰兩年第1期觀測差值; (b)-(d) 基于相鄰兩年第2期觀測差值. 洋紅色線為研究區(qū)域邊界，灰色線為省界，黑色線為陸海邊界. 單位為10-8 m/s2

圖 5 華北地區(qū)由Slepian方法確定的120階2年尺度重力場變化空間分布：(a) 基于2013.04和2011.03的觀測差值; (b) 基于2013.08和2011.08的觀測差值. 洋紅色線為研究區(qū)域邊界，灰色線為省界，黑色線為陸海邊界. 單位為10-8 m/s2

本文還利用Slepian方法計算了華北地區(qū)2011–2013年間2年尺度的GRACE重力變化（2013.04–2011.03和2013.07–2011.08），如圖6所示。

圖6華北地區(qū)由Slepian方法得到的GRACE重力變化空間分布，其中(a)、(b)分別為2013.04–2011.03和2013.07–2011.08結(jié)果. 洋紅色線為研究區(qū)域邊界，灰色線為省界，黑色線為陸海邊界. 單位為10-8 m/s2

對于2年尺度的重力變化，地面結(jié)果（圖5）與GRACE結(jié)果（圖6）表現(xiàn)出較為一致的“南減北增”空間分布特征。進(jìn)一步對比圖5和圖6，可發(fā)現(xiàn)2年尺度的地面與GRACE衛(wèi)星結(jié)果差異最大區(qū)域出現(xiàn)在石家莊–衡水連線一帶，2011–2013年區(qū)域內(nèi)地面結(jié)果為增大趨勢，GRACE結(jié)果為減少趨勢。該區(qū)域長期存在嚴(yán)重的淺層、深層地下水超量開采以及地表沉降問題，這是影響區(qū)域內(nèi)重力場長期變化的兩個主要因素。對地面測量而言，地下水虧損會令重力觀測減小，地表沉降會導(dǎo)致重力觀測增加。根據(jù)Shen等(2015)的地面觀測結(jié)果，2009–2013年在河北省南部（包含石家莊–衡水連線一帶）由地表沉降導(dǎo)致的重力增加效應(yīng)要大于地下水虧損帶來的減小效應(yīng)；由于2009–2013年間石家莊–衡水連線一帶地表沉降速度基本維持不變(約-0.05 m/a)，我們推斷2011–2013年區(qū)域內(nèi)地表沉降的影響仍然大于地下水虧損的影響，因此地面重力結(jié)果呈增大趨勢。對GRACE觀測而言，在衛(wèi)星軌道處地面沉降導(dǎo)致的重力變化已大幅衰減，此時地下水虧損導(dǎo)致的重力減小效應(yīng)占主導(dǎo)地位，因此GRACE結(jié)果整體為減少趨勢。

本文的創(chuàng)新及取得的成果如下：(1)              首次嘗試?yán)?/span>Slepian局部譜分析方法和多期地面重力觀測確定更高空間分辨率的時變重力場模型，討論給出了相應(yīng)構(gòu)建流程。(2)              提出了Slepian區(qū)域重力場建模參數(shù)的優(yōu)化方法。(3)              獲得了2011–2013年華北地區(qū)120階年尺度時變重力場模型，并驗證了模型的可靠性。

本文研究結(jié)果可為區(qū)域重力場建模提供新的參考，也可為華北地區(qū)水資源變化監(jiān)測、構(gòu)造活動分析以及地震風(fēng)險性評估等研究提供高分辨率的時變重力場模型支撐。

研究成果發(fā)表于“地球物理學(xué)報”。研究受國家自然科學(xué)基金(41974095)，中國博士后科學(xué)基金(2018M641424)，中國地震局地球物理研究所基本科研業(yè)務(wù)費專項(DQJB19B03,DQJB20X09)聯(lián)合資助。

【相關(guān)文獻(xiàn)列表】

[1].      Chen S, Zhuang J, Li X et al. 2019. Bayesian approach for network adjustment for gravity survey campaign: methodology and model test. Journal of Geodesy, 93:681-700, doi:10.1007/s00190-018-1190-7.

[2].      Shen C, Xuan S, Zou Z et al. 2015. Trends in gravity changes from 2009 to 2013 derived from ground-based gravimetry and GRACE data in North China. Geodesy and Geodynamics, 6(6):423-428, doi:10.1016/j.geog.2015.08.001.

[3].      Simons F J. 2010. Slepian functions and their use in signal estimation and spectral analysis, in: Freeden W, Nashed M Z, Sonar T (Eds.), Handbook of Geomathematics. Springer, Heidelberg, Germany, pp. 891-923.

[4].      Simons F J, Dahlen F A. 2006. Spherical Slepian functions and the polar gap in geodesy. Geophysical Journal International, 166(3):1039-1061, doi:10.1111/j.1365-246X.2006.03065.x.

[5].      Simons F J, Dahlen F A, Wieczorek M A. 2006. Spatiospectral Concentration on a Sphere. SIAM Review, 48(3):504-536, doi:10.1137/s0036144504445765.

[6].      陳石, 徐偉民, 王謙身. 2017. 應(yīng)用Slepian局部譜方法解算中國大陸重力場球諧模型. 測繪學(xué)報, 46(8):952-960, doi:10.11947/j.AGCS.2017.20150542.

[7].      韓建成, 陳石, 盧紅艷, 徐偉民. 2021. 基于Slepian 方法和地面重力觀測確定時變重力場模型：以2011–2013 年華北地區(qū)數(shù)據(jù)為例. 地球物理學(xué)報, 64(5):1542-1557, doi: 10.6038/cjg2021O0240.

【作者簡介】

韓建成2012年畢業(yè)于武漢大學(xué)測繪學(xué)院，獲固體地球物理學(xué)博士學(xué)位。2014年至2017年于臺灣交通大學(xué)從事博士后研究，2018年進(jìn)入中國地震局地球物理研究所從事博士后研究，現(xiàn)為第四研究室助理研究員。目前主要研究方向為重力場恢復(fù)、時變重力場應(yīng)用及高性能計算等。