王建強(qiáng) 孫云龍

1 東華理工大學(xué)測(cè)繪工程學(xué)院，南昌市廣蘭大道418號(hào)，330013

高程基準(zhǔn)的統(tǒng)一是大地測(cè)量的長期目標(biāo)[1]，對(duì)工業(yè)、經(jīng)濟(jì)及大型工程建設(shè)都具有重要的支撐作用。隨著GNSS技術(shù)的廣泛應(yīng)用及我國北斗系統(tǒng)的全面建成，構(gòu)建高精度高程基準(zhǔn)成為高程測(cè)量實(shí)時(shí)化建設(shè)需要突破的關(guān)鍵技術(shù)。發(fā)達(dá)國家近幾十年不斷精化本國大地水準(zhǔn)面，美國最新大地水準(zhǔn)面模型Geoid18[2]的分辨率優(yōu)于2 km，內(nèi)符合精度優(yōu)于2 cm；澳大利亞發(fā)布的AUSGeoid2020內(nèi)符合精度達(dá)到亞cm級(jí)[3]；加拿大同美國正合作構(gòu)建Geoid2022，屆時(shí)可全面覆蓋加拿大區(qū)域[4]。我國的大地水準(zhǔn)面研究技術(shù)不論是在理論還是工程應(yīng)用方面均已走在世界前列[5-6]，李建成[1]在2020-10中國測(cè)繪學(xué)會(huì)年會(huì)上指出，中國似大地水準(zhǔn)面的擬合精度已達(dá)2～3 cm，但構(gòu)建高精度(似)大地水準(zhǔn)面仍有提升空間，還需開展更多研究分析。本文研究構(gòu)建市級(jí)區(qū)域似大地水準(zhǔn)面的虛擬球諧方法，并以南昌市為例進(jìn)行計(jì)算分析。

1 位模型計(jì)算大地水準(zhǔn)面

地球外部引力位的球諧表達(dá)式[7]為：

(1)

(2)

式中，J2n為系數(shù)值，當(dāng)J2已知時(shí)，其他數(shù)值可以計(jì)算得出；P2n(θ)為勒讓德函數(shù)。將式(1)減去式(2)可得擾動(dòng)位重力異常位T：

T(r,θ,λ)=V(r,θ,λ)-U(r,θ,λ)

(3)

(4)

(5)

(6)

式中，γ為正常重力值，可通過公式直接計(jì)算獲取。

利用式(6)計(jì)算出的大地水準(zhǔn)面有2個(gè)誤差源，分別為位系數(shù)誤差和截?cái)嗟絅max階所引起的誤差，后者被稱為模型截?cái)嗾`差，可采用重力異常進(jìn)行估算。大地水準(zhǔn)面上的重力異常表達(dá)式為:

(7)

(8)

由此可得：

(9)

將式(9)代入式(6)可得：

(10)

擾動(dòng)引力分量的截?cái)嗾`差表達(dá)式為：

(11)

式中，R為地球平均半徑；Cn可利用Moritz兩分量模型和Lapp參數(shù)進(jìn)行計(jì)算：

(12)

其中，C2=7.5 mGal2，計(jì)算階次的上限為100 000，R=6 371 km，結(jié)果如圖1所示。由圖可知，要達(dá)到cm級(jí)精度的大地水準(zhǔn)面，截?cái)嚯A次要達(dá)到2 000以上。

圖1 大地水準(zhǔn)面模型截?cái)嗾`差Fig.1 Geoid model truncation error

2 虛擬球諧理論與方法

球冠諧模型是構(gòu)建區(qū)域大地水準(zhǔn)面模型的典型方法[8]。在球冠坐標(biāo)系下，球冠半徑為θ0，任一點(diǎn)的坐標(biāo)為(r,θ,λ)，θ為余緯，λ為經(jīng)度。在球冠諧分析中，地球外部重力場(chǎng)位滿足Laplace方程，同樣可用分離變量法獲得方程的解[9]。由于球冠諧分析中，余緯的邊界條件[9]為：

T(r,θ0,λ)=f(r,λ)

(13)

(14)

可以看出，兩式等號(hào)右端的函數(shù)均與θ無關(guān)，而在球諧分析中，僅有締合勒讓德函數(shù)與θ有關(guān)。式(13)和式(14)的基函數(shù)可以通過以下2個(gè)方程分別滿足：

(15)

(16)

由于其經(jīng)度范圍和球諧分析中的經(jīng)度范圍相同，因此對(duì)應(yīng)的本征值m為整數(shù)變量。當(dāng)給定式(15)和式(16)的θ0時(shí)，可單獨(dú)確定一組對(duì)應(yīng)m的n值序列，此時(shí)的n為非整數(shù)，非整階締合勒讓德函數(shù)的特征值可通過Muller方法計(jì)算得到[10]。通過球冠來確定非整階勒讓德函數(shù)序列[11]時(shí), 需要2個(gè)正交基數(shù)，令k為n值序列的下標(biāo)，并定義k-m為奇數(shù)時(shí)，采用式(15)獲取的本征值序列；k-m為偶數(shù)時(shí)采用式(16)獲取的本征值序列?？梢钥闯?，球冠諧函數(shù)描述的地球重力場(chǎng)是局部區(qū)域的，與全球區(qū)域的球諧函數(shù)存在差異：在求解關(guān)于余緯的偏微分方程時(shí)，球諧分析中方程的本征值是整數(shù)，而球冠諧分析中的本征值是非整數(shù)。擾動(dòng)位球冠諧展開式為：

(17)

由于受非整階締合勒讓德函數(shù)的限制，球冠諧模型很難達(dá)到高分辨率，球冠半徑較小(如5°以下)的函數(shù)變化更為迅速。圖2為球冠半徑為5°，m=0、1、2、4時(shí)的函數(shù)值，可以看出，能計(jì)算出的零根值數(shù)量在10個(gè)左右。Haines[13]利用逼近方法改進(jìn)算法，但仍有較大誤差。利用球冠諧模型逼近區(qū)域大地水準(zhǔn)面可達(dá)到cm級(jí)的逼近精度，但隨著區(qū)域的擴(kuò)大及地形的復(fù)雜化，球冠諧模型仍受階次擴(kuò)展的限制。為克服這一限制，借用球冠諧映射方法，將球冠坐標(biāo)系變換后采用整階次締合勒讓德函數(shù)進(jìn)行計(jì)算[14]。

圖2 締合勒讓德函數(shù)值Fig.2 Legendre values

確定逼近區(qū)域后，根據(jù)球冠諧理論設(shè)計(jì)球冠半徑和球冠中心，在保持相當(dāng)精度的前提下，將球冠諧函數(shù)改進(jìn)為虛擬球諧函數(shù)。將余緯θ的定義域(0，θ0)映射到(0，π)上，用整階次締合勒讓德函數(shù)代替非整階締合勒讓德函數(shù)。虛擬球諧方法需要將原坐標(biāo)系轉(zhuǎn)換到新坐標(biāo)系中：

(18)

3 數(shù)值實(shí)驗(yàn)

以南昌市為例進(jìn)行實(shí)驗(yàn)(圖3)，數(shù)據(jù)點(diǎn)由主要交通路線隨機(jī)采樣獲得。通過EGM2008模型計(jì)算該區(qū)域內(nèi)階次分別為2～60、61～360、361～1 000、1 001～2 100的大地水準(zhǔn)面，計(jì)算結(jié)果如圖4所示?？梢钥闯?，60階次內(nèi)的重力場(chǎng)模型計(jì)算的大地水準(zhǔn)面起伏范圍達(dá)到2 m，61～360階次范圍內(nèi)的大地水準(zhǔn)面起伏范圍達(dá)到40 cm，360階次以上的大地水準(zhǔn)面起伏范圍大約為20 cm，此時(shí)的起伏范圍已經(jīng)較小。對(duì)圖3中的采樣點(diǎn)進(jìn)行數(shù)值計(jì)算，得到剩余大地水準(zhǔn)面數(shù)據(jù)范圍為-42～-30 cm，具體見圖5。在此基礎(chǔ)上，增加1 cm 的白噪聲誤差。

圖3 數(shù)據(jù)采集范圍Fig.3 Data collection range

圖4 不同階次大地水準(zhǔn)面起伏情況Fig.4 The geoid undulates in different orders

圖5 去除360階次重力場(chǎng)模型后的大地水準(zhǔn)面數(shù)據(jù)Fig.5 The geoid after removing the 360-ordergravitational field model

通過重力場(chǎng)模型濾波后得到新的觀測(cè)值，然后將地理坐標(biāo)轉(zhuǎn)換為球冠坐標(biāo)。實(shí)驗(yàn)數(shù)據(jù)范圍大約在1°左右，其中半徑為0.5°，添加0.5°的過渡帶。根據(jù)虛擬球諧理論設(shè)置虛擬半徑[15]為1°，模型階次設(shè)定為13，利用式(18)將球冠坐標(biāo)轉(zhuǎn)換為新的虛擬坐標(biāo)，并在新坐標(biāo)系下利用球諧函數(shù)構(gòu)建模型。為突出虛擬球諧理論的優(yōu)越性，加入多項(xiàng)式擬合和神經(jīng)網(wǎng)絡(luò)模型進(jìn)行對(duì)比，結(jié)果如圖6所示。

圖6 誤差分布Fig.6 Error distribution

實(shí)驗(yàn)結(jié)果表明，虛擬球諧方法的擬合精度最高，達(dá)到1 cm，RMS值和STD值是4種方法中最小的，都只有0.308 cm；其次是神經(jīng)網(wǎng)絡(luò)模型，擬合精度基本優(yōu)于2 cm；4次多項(xiàng)式的擬合誤差較大，接近5 cm。具體實(shí)驗(yàn)結(jié)果如表1所示。

表1 精度表

4 結(jié) 語

本文研究構(gòu)建市級(jí)大地水準(zhǔn)面的虛擬球諧理論及方法，并以南昌市為例，對(duì)比分析了虛擬球諧方法、多項(xiàng)式擬合方法及神經(jīng)網(wǎng)絡(luò)模型。從實(shí)驗(yàn)結(jié)果可以看出，虛擬球諧方法的RMS值只有0.308 cm，低于多項(xiàng)式擬合和神經(jīng)網(wǎng)絡(luò)模型，表明虛擬球諧方法偏差更小、效果更好。本文僅采用重力場(chǎng)位模型數(shù)據(jù)對(duì)觀測(cè)值進(jìn)行簡(jiǎn)單的移去和恢復(fù)操作，沒有綜合利用重力數(shù)據(jù)和地形數(shù)據(jù)，因此未來仍需作進(jìn)一步探索研究。