管雅慧,王坦, 胡順強,楊振宇

(1.首都師范大學 地球空間信息科學與國際化示范學院,北京 100048;2.首都師范大學 資源環(huán)境與旅游學院,北京 100048;3.中國地震臺網(wǎng)中心,北京 100045)

0 引 言

全球定位系統(tǒng)(GPS)自建立以來累積了大量觀測數(shù)據(jù),為分析非線性構(gòu)造運動的主要來源提供可靠數(shù)據(jù),對于研究GPS時間序列的噪聲特性有著重要意義. 通過分析非構(gòu)造信息可以提高基準站速度估計的準確性,為進一步研究該地區(qū)地殼形變模式提供理論依據(jù). 針對時間序列中普遍存在有色噪聲的現(xiàn)象,Williams等[1]對9個不同時間跨度的GPS全球解及局域網(wǎng)解進行噪聲分析,判定在GPS全球解中白+閃爍噪聲(WN+FN)模型是最為合適的噪聲模型組合,GPS局域網(wǎng)解則難以判定噪聲模型的優(yōu)劣,測站相關(guān)因素對其影響較大. 丁開華等[2]利用4種噪聲模型分析川滇地區(qū)52個GPS基準站三年坐標時間序列的噪聲特性,認為白+冪律噪聲(WN+PL)模型或WN+FN噪聲模型最能體現(xiàn)該地區(qū)基準站的噪聲特性.張風霜[3]采用6種噪聲模型對云南范圍內(nèi)25個基準站四年的坐標時間序列進行研究,發(fā)現(xiàn)三個坐標分量上的主要噪聲特性都體現(xiàn)為WN+FN和WN+PL噪聲模型,兩個模型在水平方向上共占76%,在垂直方向上共占100%,分布在紅河斷裂帶兩側(cè)垂直方向上的最優(yōu)噪聲模型呈現(xiàn)一定的地域性. 另外,賀小星等[4]指出隨著時間跨度的增加,全球范圍內(nèi)基準站時間序列下的最優(yōu)噪聲模型中隨機漫步噪聲(RW)模型的比例會變高. 綜合上述研究表明: 1)不同地區(qū)GPS基準站各分量所體現(xiàn)的噪聲特性并不相同. 2)使用不同的噪聲模型評價準則和方法判定最優(yōu)模型會得到不同的結(jié)果. 3)時間序列的最優(yōu)噪聲模型會隨著時間跨度的增加發(fā)生一定的變化.若無法有效檢測出坐標分量時間序列中存在的有色噪聲,將會極大低估觀測噪聲帶來的誤差影響,無法有效分離噪聲與形變信號. 因此,正確有效判定基準站噪聲模型的優(yōu)劣是十分必要的.

1 云南地區(qū)GPS基準站數(shù)據(jù)處理

1.1 數(shù)據(jù)準備及解算策略

自“十一五”計劃期間啟動了中國陸態(tài)網(wǎng)絡(luò)工程以來,在云南地區(qū)陸續(xù)建成27個GPS基準站,如圖1所示. 此次解算使用的GPS數(shù)據(jù)是陸態(tài)網(wǎng)絡(luò)27個GPS基準站自2011-01-01到2017-10-14的連續(xù)觀測資料,聯(lián)合15個中國周邊國際衛(wèi)星全球定位服務(wù)(IGS)基準站(BJFS、DAEJ、GUAM、GUAO、HYDE、IISC、KIT3、LHAZ、NVSK、PIMO、SELE、TCMS、TWTF、TSKB、 YIBL)一同參與解算.

圖1 云南地區(qū)27個GPS基準站分布圖

本次數(shù)據(jù)處理選擇的是GAMIT/GLOBK 10.61版本數(shù)據(jù)處理軟件,是世界上應(yīng)用最為廣泛的高精度GPS后處理軟件之一.在GAMIT基線解算過程中,控制參數(shù)設(shè)置不同會給基線結(jié)果精度帶來較大影響,因此正確恰當?shù)卦O(shè)置控制參數(shù)顯得尤為重要. 表1是本次GAMIT基線解算控制參數(shù)設(shè)置表.

表1 解算參數(shù)設(shè)置表

完成GAMIT基線解算后,將基線松弛解h文件與來自IGS數(shù)據(jù)中心的全球解h文件數(shù)據(jù)進行合并,獲取集合所有數(shù)據(jù)的松弛解. 利用卡爾曼濾波對該松弛解進行網(wǎng)平差,得到基準站坐標時間序列,圖2所示為解算的YNLA站的原始坐標時間序列.

圖2 YNLA站點的坐標時間序列

1.2 解算精度分析

評價GAMIT/GLOBK軟件解算成果精度應(yīng)滿足以下兩個指標.

1) 單日解q文件中的標準均方根誤差(NRMS)值.GAMIT解算基線采用全組合解的形式,歷元模糊解計算中的殘差是NRMS值的來源,NRMS值可以直觀反映基線解在某一時段內(nèi)解算偏離基線解加權(quán)平均值的程度[5]. 如果NRMS值越低,表明基線解質(zhì)量越好,NRMS值正常情況下應(yīng)小于0.3. 若存在部分周跳在解算中未得到良好修復或某參數(shù)出現(xiàn)偏差,將會導致NRMS值高于0.5.其公式如下:

(1)

式中:N為解算使用的基線總量;Yi為第i天所求的基線解;Y為第i天單位時段內(nèi)目標基線的加權(quán)平均值;σi為第i天該基線的中誤差[6].

提取所有單日解q文件中NRMS值,將其繪制成表，如表2所示，NRMS值全部低于0.5,絕大部分NRMS值均處在在0.2附近. 說明基線精度達到網(wǎng)平差的要求,本次解算結(jié)果質(zhì)量較好.

表2 NRMS值分布

2) GPS觀測質(zhì)量和基線解算精度的另一重要衡量指標是單日解重復性,表現(xiàn)在各方向分量及邊長的基線重復性. 如果基線重復性越高,表明GPS觀測質(zhì)量越差,基線解算精度越低,反之觀測質(zhì)量越好,解算精度越高[7]. 基線重復性的計算公式如下:

(2)

(3)

通過固定誤差和比例誤差的方式表達單日解的基線重復性可以更為直觀地展示其解算精度:

RC=a+b·Lc，

(4)

式中:a表示基線重復性在各方向的固定誤差;b表示所求基線分別在各方向的比例誤差;Lc表示單日解基線在各方向分量上的長度.

表3 各方向的重復性統(tǒng)計結(jié)果

本次解算的基線重復性在南北方向上優(yōu)于1.40 mm+2.02×10-9×S,東西方向上優(yōu)于2.30 mm+1.68×10-9×S,垂直方向上優(yōu)于7.05 mm+2.44×10-9×S,基線長度重復性優(yōu)于1.69 mm+1.07×10-9×S(S表示基線長度,mm為單位),滿足我國相關(guān)GPS測量規(guī)范的要求.

1.3 坐標序列預處理

由于受到天氣變化、接收機故障或信號不暢等因素影響,GLOBK獲取的原始坐標時間序列中還會存在少量粗差.為降低粗差對后續(xù)噪聲分析結(jié)果的影響,采取四分位距(IQR)判別準則對時間序列的粗差進行剔除,IQR判別準則原理如下:

IQR=Q2-Q1.

(5)

異常值探測區(qū)間為

[Q1-1.5·IQR,Q2+1.5·IQR].

(6)

IQR準則假設(shè)該時間序列數(shù)據(jù)滿足標準正態(tài)分布形式,按照從小到大的順序依次排列.式中:Q1表示為最靠近1/4處的下分位值;Q2表示為最靠近3/4處的上分位值. 若出現(xiàn)IQR統(tǒng)計值偏離3倍以上坐標分量殘差的情況,則剔除該單日解坐標殘差[9].

由于原始坐標時間序列進行過預處理,坐標時間序列中的數(shù)據(jù)缺失是無法避免的. 因此需要將序列中缺失間斷的部分插值擬合. 間隔時間不同的缺失數(shù)據(jù)需要不同的插值方法,三次樣條法適用于間隔時間小于3天的缺失時間序列,線性插值法適用于間隔時間大于3天的缺失時間序列,避免坐標時間序列原有的變化趨勢被破壞[5],圖3為預處理后的YNLA站的坐標時間序列圖.

圖3 插值擬合后的YNLA站點的坐標時間序列

2 GPS時間序列的噪聲分析

2.1 GPS時間序列函數(shù)模型

為了分析GPS基準站坐標時間序列中各種運動參數(shù)信息,對各基準站運動模型建立如下公式[10]:

y(ti)=a+bti+csin(2πti)+dcos(2πti)+

esin(4πti)+fcos(4πti)+

H(ti-Tkj)+vi

(7)

2.2 最優(yōu)有色噪聲模型評價準則

以往評估噪聲模型的傳統(tǒng)方法多以最大似然估計法(MLE)和功率譜分析法為主[11]. 功率譜法可以分析經(jīng)過擬合處理后的殘差時間序列噪聲模型,通過譜指數(shù)所在區(qū)間判定主要的噪聲類型,但該方法對低頻噪聲分量估計不夠準確,計算結(jié)果存在較大誤差. 通過比較多種噪聲模型組合分析得到的MLE值,選取得到最大MLE值的噪聲模型組合為最優(yōu)有色噪聲模型. 但MLE法對噪聲模型未知參數(shù)無法起到較好的約束作用,當參與噪聲分析的未知參數(shù)越多時,MLE值越大,因此單純對比MLE值并不能準確反映出坐標時間序列的噪聲特性. 針對MLE方法的缺陷,Langbein[12]認為設(shè)定經(jīng)驗值做零假設(shè)分析可以有效區(qū)分不同的噪聲模型組合優(yōu)劣程度,得出最終的最優(yōu)噪聲模型. 然而經(jīng)驗值設(shè)定的不確定性在較大程度上會影響結(jié)果,需要進行更多的研究. 通過分析以上模型評價方法的優(yōu)缺點,本文評價最優(yōu)噪聲模型的方法為赤池信息量準則(AIC)和貝葉斯信息準則(BIC)[13-14]. 根據(jù)AIC/BIC值判定不同噪聲模型的優(yōu)劣來確定最優(yōu)模型,其基本原理如下:

AIC=-2ln(L)+2K,

(8)

BIC=-2ln(L)+Kln(n) .

(9)

式中:L為目標模型中的似然函數(shù)；K為該模型內(nèi)參數(shù)個數(shù)；n為觀測值總量.通過式(8)和(9)計算可知,L值越大,對應(yīng)的AIC/BIC值越小,說明對應(yīng)的模型與時間序列越相近[15]. 在大部分情況下AIC/BIC最小值對應(yīng)的噪聲模型基本一致,在少數(shù)特殊情況下,AIC/BIC最小值對應(yīng)的噪聲模型并不完全相同.由于BIC值是依據(jù)一系列參數(shù)模型擬合后的真實模型,因此應(yīng)選取BIC值較小的模型組合作為最優(yōu)噪聲模型.

3 時間序列跨度分析實例

3.1 時間跨度對最優(yōu)噪聲模型的影響

本文以27個基準站時間跨度分別為3 a、5 a、7 a的GPS時間序列為研究對象,采用六種有色噪聲模型組合做噪聲估計,白噪聲(WN)、WN+FN、PL、閃爍+隨機漫步噪聲(FN+RW)、白+閃爍+隨機漫步噪聲(WN+FN+RW)、廣義高斯馬爾科夫噪聲模型 (GGM),比較各基準站在不同跨度不同噪聲模型下的AIC/BIC值,獲取各站在N、E、U三個方向的最優(yōu)噪聲模型. 圖4為不同時間跨度下N、E、U各方向最優(yōu)噪聲模型分布圖.

圖4 不同時間跨度下N、E、U各方向最優(yōu)噪聲 模型分布圖

由圖4可知:根據(jù)不同時間跨度各分量上最優(yōu)噪聲模型的對比,時間序列三個方向上的最佳噪聲模型都會經(jīng)歷從發(fā)散到整體一致性逐漸提高的過程. 當時間跨度為3 a時,在N、E、U方向上表現(xiàn)的噪聲特性都是以PL模型為主(分別約占51.9%、48.2%、85.2%),WN+FN模型為輔(分別約占40.7%、44.4%、22.2%). 少數(shù)測站在水平方向上的主要噪聲特性表現(xiàn)為GGM模型. 當時間序列跨度積累到5 a時,E方向上最優(yōu)噪聲模型所占比例變化趨勢并不明顯,其趨勢主要表現(xiàn)在N方向和U方向上,PL模型所占比例有所下降(分別約占25.9%、66.7%),WN+FN模型所占比例有所上升(分別約占59.3%、22.2%). 三個分量上最優(yōu)噪聲模型的類型減少,逐漸趨于一致. 當時間跨度為7 a時,E方向上最優(yōu)噪聲模型比例保持穩(wěn)定,N方向和U方向上PL模型比例持續(xù)下降(分別約占37.0%、77.8%),WN+FN模型比例持續(xù)上升(分別約占70.4%、33.3%). 同時水平方向上主要噪聲特性表現(xiàn)為WN+FN+RW模型的基準站個數(shù)也有所增加,若要進一步分析該地區(qū)基準站的隨機漫步噪聲(RW)特性還需要積累更長的時間序列.

3.2 最優(yōu)噪聲模型在不同時間跨度下的演化過程

為了探究時間跨度對建立噪聲模型的影響,不同時間跨度下各個基準站N方向的最優(yōu)噪聲模型演化過程如表4所示.

表4 不同時間跨度下噪聲模型的演化過程(以N方向為例)

當時間跨度由3 a增大到5 a時,最優(yōu)噪聲模型發(fā)生改變的基準站約占37%,說明部分跨度為3 a的時間序列表現(xiàn)的噪聲特性并不穩(wěn)定,噪聲中長周期分量還未準確探測,結(jié)果具有發(fā)散性;時間跨度從5 a到7 a的過程中約11.1%的基準站最優(yōu)噪聲模型發(fā)生改變,說明大部分跨度為5 a的時間序列最優(yōu)噪聲模型并未隨著時間跨度的改變而改變,結(jié)果已具備一定的可靠性. 針對E方向最優(yōu)噪聲模型演化結(jié)果分析可知,約18.5%的基準站在時間跨度增大到5 a后最優(yōu)噪聲模型發(fā)生改變;當時間跨度增大到7 a時約88.9%的基準站最優(yōu)噪聲模型保持不變. U向最優(yōu)噪聲模型統(tǒng)計結(jié)果表明,分別有25.9%和14.8%的基準站在5 a和7 a時間跨度后最優(yōu)噪聲模型發(fā)生變化. 通過對三個方向噪聲模型在不同時間跨度的演化規(guī)律研究可知:當時間跨度大于5 a時,基準站噪聲模型逐漸趨于穩(wěn)定. 時間跨度積累越長,更有利于未知噪聲的探測,噪聲模型估計的準確度越高.因此要想獲得較為可靠的時間序列噪聲特性信息,至少需要積累超過5 a的GPS時間序列.

4 結(jié) 論

本文基于云南地區(qū)基準站長達7 a的GPS有效觀測數(shù)據(jù),求出27個GPS基準站在各方向上的最優(yōu)有色噪聲模型,通過對不同時間跨度3 a、5 a、7 a的GPS時間序列噪聲模型演化規(guī)律分析,結(jié)果表明:本次解算的基線重復性在水平方向上優(yōu)于2.30 mm+1.68×10-9×S,垂直方向上優(yōu)于7.05 mm+2.44×10-9×S,單日解的標準化均方差均處在0.2左右.時間序列噪聲模型建立會受到時間跨度的影響,增大時間跨度會引起N方向和U方向上最優(yōu)噪聲所占模型比例的改變,主要體現(xiàn)在PL模型所占比例有所下降,WN+FN模型所占比例有所上升. 當時間跨度大于5 a時,基準站噪聲模型的穩(wěn)定性顯著提高,同時增加未知噪聲分量出現(xiàn)的可能性.