王文超 王夏莉 王 斌 耿曉燕 劉煥玲
1 國家自然資源部大地測量數(shù)據(jù)處理中心,西安市友誼東路334號,7100542 中國測繪科學(xué)研究院,北京市蓮花池西路28號,100036
帶有齒輪系統(tǒng)及量測螺旋系統(tǒng)的LCR-G型和Burris型相對重力儀由于其特殊的結(jié)構(gòu),會存在周期誤差,該誤差主要由齒輪偏心、齒輪齒距分度誤差、量測螺桿螺紋距誤差等機制產(chǎn)生[1]。短周期項對重力基準(zhǔn)網(wǎng)解算精度的影響平均約在3 μGal以下,而中長周期項的影響最大超過20 μGal[2-3]。因此,重力基準(zhǔn)網(wǎng)的解算必須考慮周期項的影響。目前,重力基準(zhǔn)網(wǎng)對周期項的解算主要是利用其周期項的展開式。一般情況下需考慮7個周期項,但是當(dāng)7個周期項全部考慮時,由于參數(shù)過多會發(fā)生過擬合現(xiàn)象,且7個周期項全部解算精度要達(dá)到2 μGal以上大約需要200多個觀測值[3-4],但在實際觀測中很難保證每臺儀器均獲得200個以上的觀測量。因此,相對重力儀周期項的選取是重力基準(zhǔn)網(wǎng)解算的一個難題。
目前相對重力儀周期項的選取主要采用以下2種方法:1)根據(jù)誤差傳播定律,得出周期項的取舍標(biāo)準(zhǔn)[5],然后對周期項進(jìn)行篩選;2)采用逐步回歸分析法[6],利用多元回歸方程中各變量方差的貢獻(xiàn)進(jìn)行顯著性檢驗,根據(jù)檢驗結(jié)果對周期項進(jìn)行保留或剔除。但這2種方法在篩選過程中均需要大量的觀測量來支撐。在機器學(xué)習(xí)領(lǐng)域,有學(xué)者使用正則化方法[7-8]對多項式回歸方程的冗余項進(jìn)行剔除,結(jié)果表現(xiàn)良好,能有效防止過擬合現(xiàn)象。因此,本文基于2000國家重力基本網(wǎng)及2020最新國家重力基準(zhǔn)網(wǎng)的觀測數(shù)據(jù),使用正則化方法,以LCR-G型相對重力儀為例,研究國家重力基準(zhǔn)網(wǎng)中相對重力儀周期項參數(shù)的確定方法。
相對重力儀周期項的函數(shù)模型[9]為:
(1)
式中,An為各周期項的振幅,φn為各周期項的初相,R為儀器讀數(shù),Tn為相對重力儀的傳動周期,一般情況下,n≤7。
由于該函數(shù)模型的形式不適合平差過程中的誤差方程組建,因此,對其進(jìn)行推導(dǎo)變換。令Xn=Ancosφn、Yn=-Ansinφn,將Xn和Yn作為誤差參數(shù)代入式(1)得:
(2)
在重力基準(zhǔn)網(wǎng)平差解算過程中,相對重力儀的參數(shù)結(jié)果分為以下幾種情況:
1)如果Xn=0、Yn=0,則振幅An為0,該周期項可以忽略;
重力基準(zhǔn)網(wǎng)平差采用間接平差模型,將經(jīng)過固體潮改正、海潮改正、氣壓改正、儀器高改正、零漂改正后測段起、止點的最后觀測值和儀器讀數(shù)作為觀測量,其誤差方程為[4]:
(3)
式中,gi、gj分別為測站i、j點平差后的重力值,gRZi、gRZj分別為測站i、j點經(jīng)過5項改正的相對聯(lián)測最后觀測值,Ri、Rj分別為儀器在測站i、j點的觀測讀數(shù),CK為重力儀的M次多項式格值函數(shù)的K次格值改正因子,Xn、Yn為重力儀周期誤差參數(shù),Tn為周期誤差的周期,其周期值大小見表1,1個計數(shù)單位(1格)相當(dāng)于mGal。
表1 LCR-G型重力儀傳動周期
正則化函數(shù)模型[7-8]公式為:
Rλ(β)=
(4)
式中,N為觀測量的個數(shù),β為待求量,xi和yi為觀測量,λ為懲罰因子,α為正則化參數(shù),Pα(β)為懲罰項,其表示形式如下:
(5)
式中,p為待求參數(shù)的個數(shù)。α=1時,式(5)為Lasso算法模型;α=0時,式(5)為Ridge算法模型;α∈(0,1)時,式(5)為ElasticNet算法模型。
通過調(diào)整懲罰因子λ和正則化參數(shù)α,確定出剔除周期項的最優(yōu)參數(shù)值并對無效的周期項進(jìn)行剔除;然后再利用周期項取舍標(biāo)準(zhǔn),剔除誤差超限的周期項,從而完成對周期項參數(shù)的篩選;最后通過平差對周期項進(jìn)行精確解算,計算出儀器周期項的振幅及初相,并對其結(jié)果進(jìn)行分析。
根據(jù)誤差傳播定律,周期項的取舍標(biāo)準(zhǔn)為[2-3]:
(6)
式中,An為周期項的振幅,MA和Mφ分別為周期項的振幅中誤差和相位中誤差[4]:
(7)
式中,MX、MY分別為周期項的點位中誤差,A為對應(yīng)周期的振幅。
當(dāng)初相φ的求定誤差大于29°或者振幅誤差比小于2時,考慮剔除相應(yīng)的周期項。
本文以2000國家重力基本網(wǎng)和2020最新國家重力基準(zhǔn)網(wǎng)中G796、G922兩臺LCR-G型相對重力儀觀測數(shù)據(jù)為例,對其周期項進(jìn)行參數(shù)解算,相應(yīng)的觀測量統(tǒng)計信息見表2。由表2可以看出,兩臺儀器的觀測量個數(shù)均比待求量個數(shù)多,保證了足夠的多余觀測量來解算相對重力儀的周期參數(shù);最大段差、最小段差、平均段差分別在同一量級,可以保證所選觀測數(shù)據(jù)具有同一性。
表2 重力儀觀測量統(tǒng)計
首先,不考慮周期項對重力基準(zhǔn)網(wǎng)進(jìn)行平差,將得到的各個重力點上的重力值及所有儀器的格值一次項作為輸入數(shù)據(jù),組建由G796、G922兩臺儀器周期項構(gòu)成的多元線性方程組,并采用正則化算法進(jìn)行解算,正則化參數(shù)α分別取0.1、0.5、1.0,懲罰因子λ取0~100,由解算結(jié)果計算的均方誤差如圖1所示。從圖1可以看出,當(dāng)懲罰因子λ大于80時,均方誤差均開始增大;正則化參數(shù)α的取值對周期項的均方誤差幾乎沒有影響。由此得出,在解算儀器周期項時,懲罰因子λ取80,正則化參數(shù)α的值可以任意選取。
圖1 重力儀周期項均方誤差Fig.1 The mean square error of periodicterms of gravity meters
采用正則化方法,懲罰因子λ取80,正則化參數(shù)α分別取0.1、0.5、1.0,解算的G796、G922儀器周期參數(shù)結(jié)果見圖2??梢钥闯?,正則化方法對周期項的剔除效果明顯,正則化參數(shù)α的選取對周期項的篩選幾乎不產(chǎn)生影響。
圖2 重力儀周期項振幅Fig.2 The amplitude values of periodic terms of gravity meters
采用式(6)對由正則化方法解算的周期項再次進(jìn)行判定,去掉誤差超限的周期項,最終的篩選結(jié)果見表3。
表3 重力儀周期項參數(shù)篩選結(jié)果
由于固定線性項和不固定線性項的標(biāo)定結(jié)果基本一致[10],因此,采用不固定線性項的方式來求解標(biāo)定參數(shù)。將表3中篩選出的有效周期項確定的參數(shù)代入國家重力基準(zhǔn)網(wǎng)的誤差方程,重新進(jìn)行平差得到周期項的振幅和相位,結(jié)果見表4。由表4可知,通過正則化方法確定的周期項,振幅中誤差優(yōu)于4.2 μGal,相位中誤差優(yōu)于27.7,解算結(jié)果精度均在周期項取舍標(biāo)準(zhǔn)要求以內(nèi)。2000國家重力基本網(wǎng)G796重力儀周期項振幅最大可達(dá)12 μGal,G922重力儀周期項振幅最大可達(dá)8.0 μGal;2020最新國家重力基準(zhǔn)網(wǎng)G796重力儀周期項振幅最大可達(dá)10.7 μGal,G922重力儀周期項振幅最大可達(dá)7.5 μGal。兩臺儀器的振幅在近20 a內(nèi)均有所下降,可能是由于儀器出廠時間較長、部件老化引起的;相同儀器在近20 a內(nèi)所確定的周期并不相同,分析認(rèn)為是相對重力儀解算的周期項的振幅和相位大小與采用的觀測數(shù)據(jù)及儀器觀測時所處的環(huán)境有關(guān)。
表4 重力儀周期項計算結(jié)果
考慮周期項前后重力基準(zhǔn)網(wǎng)平差重力值的變化見表5??梢钥闯?,考慮周期項后,對2000國家重力基本網(wǎng)的平差重力值最大影響可達(dá)6.22 μGal,對2020最新國家重力基準(zhǔn)網(wǎng)的影響可達(dá)4.73 μGal。由此可見,周期項對重力基準(zhǔn)網(wǎng)影響的數(shù)量級達(dá)到μGal,不可忽視。
表5 考慮周期項前后重力值變化的統(tǒng)計
本文基于國家重力基準(zhǔn)網(wǎng)觀測數(shù)據(jù),采用正則化方法及周期項取舍標(biāo)準(zhǔn)對G796、G922兩臺相對重力儀的周期項進(jìn)行篩選,并通過平差對周期項進(jìn)行精確解算,得出結(jié)論如下:1)使用正則化方法結(jié)合周期項取舍標(biāo)準(zhǔn),可以有效篩選出帶有齒輪系統(tǒng)及量測螺旋系統(tǒng)的相對重力儀的周期項,并且通過實驗確定,正則化方法懲罰因子λ取80最優(yōu),正則化參數(shù)α的選取對周期項幾乎沒有影響;2)周期誤差對2000國家重力基本網(wǎng)和2020最新國家重力基準(zhǔn)網(wǎng)的重力值平差結(jié)果影響最大分別可達(dá)6.22 μGal和4.73 μGal,由此可見,周期項對重力基準(zhǔn)網(wǎng)影響的數(shù)量級達(dá)到μGal;3)相對重力儀的周期會隨著觀測環(huán)境或者時間變化,最終解算結(jié)果不一致,這是由于周期函數(shù)模型與實際周期變化有差異,模型只能近似而不能真正模擬周期的變化。