聶建亮 張雪萍 郭鑫偉 王莉莉 趙文普
1 陜西測繪地理信息局,西安市友誼東路334號,710054 2 自然資源部大地測量數(shù)據(jù)處理中心,西安市友誼東路334號,710054 3 自然資源部第二地形測量隊,西安市測繪路6號,710054
為滿足社會發(fā)展需要,我國建立了新一代中國似大地水準(zhǔn)面模型CQG2000。21世紀(jì)初,陳俊勇等[1-3]和晁定波[4]基于地面重力、衛(wèi)星測高、船測重力等數(shù)據(jù),采用EGM96全球重力場,建立重力似大地水準(zhǔn)面模型。利用高程異??刂凭W(wǎng)對其進(jìn)行擬合糾正,獲取15′×15′分辨率的似大地水準(zhǔn)面模型,其中102°E以東地區(qū)中誤差小于±0.3 m,102°E以西、36°N以北地區(qū)中誤差小于±0.4 m,102°E以西、36°N以南地區(qū)中誤差小于±0.6 m。近年來,各省市陸續(xù)建立了高精度、高分辨率的區(qū)域似大地水準(zhǔn)面模型[5-11]。雖然高精度似大地水準(zhǔn)面模型能夠滿足精度要求,但高程轉(zhuǎn)換計算量巨大。以1∶10 000圖幅為例,單個圖幅的高程基準(zhǔn)轉(zhuǎn)換通常需要幾十萬到上百萬個要素點(diǎn),然而圖幅高程基準(zhǔn)轉(zhuǎn)換所需的高程精度卻較低,在遙感行業(yè),精度一般為dm級甚至m級[12-13]。區(qū)域似大地水準(zhǔn)面在地球表面空間域內(nèi)變化較緩,但在地形劇烈變化區(qū)域或地類變化過渡區(qū)域的變化較為明顯。梯度作為數(shù)值變化的表征量,可以反映數(shù)值在一定區(qū)域內(nèi)變化的大小、方向。CQG2000的梯度用格網(wǎng)點(diǎn)在經(jīng)、緯度方向上的高程異??臻g變化率表示,該數(shù)值可以量化不同比例尺圖幅內(nèi)高程異常的變化幅度,可用于分析圖幅內(nèi)平均高程異常對不同比例尺高程精度的適用性。此外,CQG2000的梯度也能反映我國大陸范圍內(nèi)重力場的空間變化特性,有助于解譯地球內(nèi)部質(zhì)量分布等相關(guān)問題。
基于此,本文利用梯度定量描繪CQG2000在全國范圍空間域內(nèi)的變化趨勢,分析CQG2000在空間域中的相關(guān)性,為建立實(shí)用化區(qū)域似大地水準(zhǔn)面模型提供基礎(chǔ)數(shù)據(jù),同時也為深入解釋似大地水準(zhǔn)面變化機(jī)理提供科學(xué)依據(jù)。
區(qū)域似大地水準(zhǔn)面受地球質(zhì)量分布的影響,在平原、丘陵及高山區(qū)的變化趨勢不同,變化劇烈區(qū)域的結(jié)果可能未達(dá)到行業(yè)應(yīng)用的精度標(biāo)準(zhǔn)。為更好地表述區(qū)域似大地水準(zhǔn)面格網(wǎng)模型變化,引入梯度函數(shù),即格網(wǎng)點(diǎn)(B,L)處的梯度:
(1)
式中,f(B,L)為區(qū)域似大地水準(zhǔn)面模型,B、L為緯度和經(jīng)度,ΔB、ΔL為緯度、經(jīng)度方向上的分辨率,m、n為格網(wǎng)模型的行、列數(shù),i、j為梯度矢量的分量方向。
對比分析多種方法發(fā)現(xiàn),三階反距離平方權(quán)差分方法計算效果較好[14]。因此,本文采用坡度的三階反距離權(quán)差分方法計算CQG2000格網(wǎng)點(diǎn)的梯度數(shù)值。
采用三階反距離平方權(quán)差分方法選取3×3的移動窗口,計算似大地水準(zhǔn)面規(guī)則格網(wǎng)點(diǎn)的梯度(圖1),即
(2)
(3)
式中,fB、fL為緯度、經(jīng)度方向上的梯度分量,ΔB、ΔL為緯度、經(jīng)度方向上的分辨率,ζ1~ζ9分別為9個格網(wǎng)點(diǎn)的高程異常值。
圖1 梯度計算示意圖Fig.1 Schematic diagram of gradient calculation
梯度方位角α為:
(4)
假設(shè)在3×3的窗口內(nèi),區(qū)域似大地水準(zhǔn)面可以表示為曲面ζk=f(Bk,Lk),且三階連續(xù)可微,按泰勒級數(shù)展開并取至三次項(xiàng),即
f(B+ΔB,L)=f(B,L)+fB(B,L)ΔB+
(5)
f(B-ΔB,L)=f(B,L)-fB(B,L)ΔB+
(6)
在8個格網(wǎng)點(diǎn)處按泰勒級數(shù)展開,整理后的fB為:
(f?B(ζB,L-ΔL)-f?B(γB,L-ΔL))+
(f?B(ζB,L+ΔL)-f?B(γB,L+ΔL))]}
(7)
令MB、ML為f?B、f?L關(guān)于B、L的三階導(dǎo)數(shù)最大范圍,其誤差具有隨機(jī)性,則有:
(8)
式中,fB一般取式(8)第1項(xiàng),第2、3項(xiàng)為數(shù)據(jù)誤差及其他誤差。
同理可得fL為:
(9)
假設(shè)區(qū)域似大地水準(zhǔn)面格網(wǎng)點(diǎn)中誤差為σζ,MB、ML中誤差相等且均為M,則σfB為:
(10)
同理可得σfL為:
(11)
由于M數(shù)值較小,因此式(10)、(11)中的第2項(xiàng)可以忽略。
根據(jù)誤差傳播定律,由式(4)、(10)、(11)推導(dǎo)出梯度方位角誤差σα為:
(12)
為更加精細(xì)地描繪CQG2000的空間域變化,在計算梯度時,對原有15′×15′分辨率的模型進(jìn)行雙線性內(nèi)插處理成5′×5′分辨率的格網(wǎng)模型(忽略內(nèi)插誤差影響)。按照§1.3中公式計算格網(wǎng)點(diǎn)梯度數(shù)值及方向,獲取全國似大地水準(zhǔn)面的梯度,如圖2所示。同時,根據(jù)表1、2中的全國梯度信息,按照CQG2000精度劃分的3個區(qū)域范圍統(tǒng)計102°E以東(區(qū)域1),102°E以西、36°N以北(區(qū)域2)和102°E以西、36 °N以南(區(qū)域3)3個區(qū)域的梯度分量精度和梯度方位角平均精度,結(jié)果見表3、圖2。
表1 CQG2000的梯度數(shù)值
表2 CQG2000的梯度方位角
表3 CQG2000的梯度精度
圖2 CQG2000的梯度分布Fig.2 The distribution of gradient for CQG2000
綜合分析表1~3和圖2可知:
1)102°E以東地區(qū)CQG2000的梯度方位角方向主要為東,量級小于0.2 m/(′),且變化平緩。其中,東北地區(qū)、華南地區(qū)梯度方向主要為東南,角度約為東向南45°;30°~40°N區(qū)間內(nèi)的梯度方位角方向?yàn)檎龞|。102°E以西地區(qū)CQG2000梯度在不同地類區(qū)域的方向不同,且梯度數(shù)值變化幅度較大。
2)我國東部、西北、西南地區(qū)的梯度分量精度σfB、σfL隨似大地水準(zhǔn)面格網(wǎng)點(diǎn)精度σζ的增大而增大,隨分辨率的增大而減小。西南地區(qū)σα誤差最大,為50.1°;西北地區(qū)σα誤差最小,為16.4°。造成上述差異的主要原因?yàn)樘荻确轿唤蔷圈姚林饕芴荻确至考捌渚鹊挠绊懀^小的梯度數(shù)值變化會導(dǎo)致較大的梯度方位角精度變化。進(jìn)一步說明,相比于梯度分量精度σfB、σfL,梯度方位角誤差σα對梯度方位角α更加敏感。
3)新疆北部與南部、青海、成都平原等地存在多個梯度擴(kuò)散中心,區(qū)域似大地水準(zhǔn)面格網(wǎng)點(diǎn)高程異常值由中心向四周逐漸增大;喜馬拉雅山脈、新疆中部地區(qū)出現(xiàn)區(qū)域似大地水準(zhǔn)面梯度會聚現(xiàn)象;河西走廊也出現(xiàn)區(qū)域似大地水準(zhǔn)面梯度會聚現(xiàn)象,甘南西部梯度方向?yàn)檎?,同時順時針旋轉(zhuǎn);寧夏南部、甘肅東北部地區(qū)區(qū)域似大地水準(zhǔn)面梯度方向?yàn)檎龞|,同時逆時針旋轉(zhuǎn)。喜馬拉雅山脈、新疆南部昆侖山脈等地的區(qū)域似大地水準(zhǔn)面梯度數(shù)值較大,最大值位于西藏墨脫縣(0.625 m/(′));青藏高原絕大部分區(qū)域的區(qū)域似大地水準(zhǔn)面梯度較小,變化較平緩。
4)CQG2000梯度變化劇烈區(qū)域與我國大山脈走向一致,說明區(qū)域似大地水準(zhǔn)面與地球內(nèi)部質(zhì)量分布密切相關(guān)。在喜馬拉雅山、昆侖山、天山、秦嶺等山脈區(qū)域,區(qū)域似大地水準(zhǔn)面變化速度較大,梯度數(shù)值變化為小-大-小。
5)CQG2000自身精度為dm級,實(shí)際分辨率為15′×15′。雖然區(qū)域似大地水準(zhǔn)面梯度精度不高,但能整體反映我國范圍內(nèi)似大地水準(zhǔn)面空間域的變化特征。此外,CQG2000梯度與全國范圍內(nèi)的垂線偏差分布[15]一致,僅在西藏東南部地區(qū)存在一定差異,可能是該區(qū)域缺乏重力等基礎(chǔ)資料、似大地水準(zhǔn)面精度較低所致。
6)由于CQG2000是利用EGM96參考重力場模型獲取高精度重力似大地水準(zhǔn)面,再通過GNSS水準(zhǔn)點(diǎn)糾正融合得到,因此重力似大地水準(zhǔn)面梯度與似大地水準(zhǔn)面梯度的總體趨勢一致,二者系統(tǒng)偏差較小,梯度大小與方向差異也較小。CQG2000梯度與高精度EGM2008重力場模型梯度在我國中東地區(qū)的差異較小,在喜馬拉雅山、昆侖山、天山、秦嶺等山脈區(qū)域存在一定差異,最大差值位于喜馬拉雅山地區(qū),約為0.18 m/(′),可能是CQG2000在該區(qū)域缺乏基礎(chǔ)資料所致。
梯度信息可反映我國似大地水準(zhǔn)面CQG2000在空間域內(nèi)的整體變化情況,描述不同范圍內(nèi)高程異常變化的量級和方向,可為建立行業(yè)所需的實(shí)用化模型提供基礎(chǔ)數(shù)據(jù)。整體上看,東北地區(qū)CQG2000梯度方位角方向?yàn)檎龞|且運(yùn)動規(guī)律與地殼水平運(yùn)動趨勢一致;西部地區(qū)梯度方位角方向隨地形分布逐漸變化,梯度劇烈變化區(qū)域與昆侖山等大山脈走向一致。本文可為進(jìn)一步開展地球動力學(xué)機(jī)理研究提供科學(xué)依據(jù)。此外,高分七號等遙感衛(wèi)星能夠開展大范圍1∶10 000立體測圖工作,滿足各行業(yè)對高精度產(chǎn)品的需求,但如何綜合利用區(qū)域似大地水準(zhǔn)面及梯度信息快速、高效地實(shí)現(xiàn)衛(wèi)星遙感影像數(shù)據(jù)處理的高程自動化轉(zhuǎn)換,是下一步的研究重點(diǎn)。