劉俊清，肖輝峰，張晨俠，丁　廣

(1. 吉林大學(xué)地球探測(cè)科學(xué)與技術(shù)學(xué)院，吉林 長(zhǎng)春 130026； 2. 吉林省地震局，吉林 長(zhǎng)春 130117； 3. 北京麥格天寶科技股份有限公司，北京 100029)

?

JLCORS網(wǎng)絡(luò)參考站坐標(biāo)時(shí)間序列噪聲類型分析

劉俊清1,2，肖輝峰3，張晨俠2，丁廣2

(1. 吉林大學(xué)地球探測(cè)科學(xué)與技術(shù)學(xué)院，吉林 長(zhǎng)春 130026； 2. 吉林省地震局，吉林 長(zhǎng)春 130117； 3. 北京麥格天寶科技股份有限公司，北京 100029)

利用吉林省連續(xù)運(yùn)行參考站網(wǎng)絡(luò)(JLCORS)近3年的觀測(cè)資料，使用GAMIT/GLOBK軟件進(jìn)行后處理，解算時(shí)用吉林省周邊IGS及陸態(tài)網(wǎng)絡(luò)站點(diǎn)作為控制點(diǎn)，獲得41個(gè)站單日松弛解的坐標(biāo)時(shí)間序列。用最大似然法擬合經(jīng)典的噪聲模型，根據(jù)噪聲模型是否能夠取得似然函數(shù)的最大值，獲得不同噪聲模型及其組合的參數(shù)估計(jì)。結(jié)果表明，JLCORS全部參考站坐標(biāo)時(shí)間序列噪聲含有冪律噪聲及白噪聲，可用二者的組合噪聲模型來(lái)描述，不含帶通濾波噪聲及各階高斯-馬爾科夫噪聲。冪律噪聲模型部分的譜指數(shù)集中分布在-2附近，表明冪律噪聲的特例隨機(jī)游走噪聲占有主導(dǎo)地位。

CORS；IGS；陸態(tài)網(wǎng)；噪聲模型

地球動(dòng)力學(xué)是建立在觀測(cè)試驗(yàn)基礎(chǔ)上的一門學(xué)科，觀測(cè)數(shù)據(jù)不可避免地存在誤差、不穩(wěn)定和不確定等非地殼構(gòu)造因素的影響[1-3]，誤差性質(zhì)及分離方法的研究一直是地學(xué)工作者的一項(xiàng)重要工作。由于地殼構(gòu)造活動(dòng)是一個(gè)非常緩慢的過(guò)程，地殼監(jiān)測(cè)時(shí)間序列跨度很長(zhǎng)，通常是幾年甚至幾十年，這就意味著觀測(cè)值中含有大量的與時(shí)間相關(guān)的誤差，而且誤差源也可能隨著時(shí)間變化[4-6]。如在多年的觀測(cè)過(guò)程中，儀器老化、更換及觀測(cè)墩標(biāo)不穩(wěn)定等，均與時(shí)間具有相關(guān)性，都會(huì)在觀測(cè)數(shù)據(jù)時(shí)間序列中產(chǎn)生誤差。隨著研究的深入，逐漸發(fā)現(xiàn)這類與時(shí)間相關(guān)的誤差的一些統(tǒng)計(jì)規(guī)律，可用白噪聲(white noise，WN)和冪率噪聲(power law noise，PLN)的組合模型進(jìn)行很好的模擬[4]。白噪聲可以通過(guò)增加觀測(cè)次數(shù)來(lái)降低影響，冪率噪聲卻無(wú)法消除,只能通過(guò)噪聲相對(duì)準(zhǔn)確的統(tǒng)計(jì)關(guān)系，在數(shù)據(jù)處理時(shí)建立數(shù)學(xué)模型來(lái)消除。

目前，很多國(guó)家建立了國(guó)家級(jí)、區(qū)域級(jí)的CORS網(wǎng)絡(luò)，主要用于測(cè)繪工作中。建站標(biāo)準(zhǔn)及其選址要求都低于地殼運(yùn)動(dòng)研究的標(biāo)準(zhǔn)，大部分站使用土層標(biāo)墩建設(shè)規(guī)范，甚至建在了建筑物頂部。這樣的觀測(cè)環(huán)境會(huì)給GNSS觀測(cè)帶來(lái)很多誤差信號(hào)，對(duì)不同用途的觀測(cè)要求，需要相應(yīng)地選擇計(jì)算模型消除誤差，這就需要對(duì)誤差性質(zhì)進(jìn)行研究，獲得其統(tǒng)計(jì)表達(dá)式，然后引入數(shù)據(jù)處理軟件。在信號(hào)研究中，誤差信號(hào)采用噪聲模型描述，本文應(yīng)用最大似然法對(duì)吉林省連續(xù)運(yùn)行參考站網(wǎng)絡(luò)(Jilin Continuously Operating Reference Station，JLCORS)參考站坐標(biāo)時(shí)間序列噪聲性質(zhì)進(jìn)行研究。

一、數(shù)據(jù)資料及噪聲分析方法

1. JLCORS概況

CORS系統(tǒng)理論源于20世紀(jì)80年代，經(jīng)過(guò)多年的發(fā)展，現(xiàn)在連續(xù)運(yùn)行參考站系統(tǒng)(CORS)能夠常年連續(xù)不斷地運(yùn)行，服務(wù)于高精度中短期天氣狀況的數(shù)值預(yù)報(bào)、變形監(jiān)測(cè)、地震監(jiān)測(cè)、地球動(dòng)力學(xué)等[7]。它不僅滿足了各種測(cè)繪、基準(zhǔn)的需求，還滿足了地震地殼形變監(jiān)測(cè)需求。目前，很多國(guó)家建立了國(guó)家級(jí)、區(qū)域級(jí)的CORS網(wǎng)絡(luò)[8]，JLCORS由吉林省測(cè)繪地理信息局建設(shè)(如圖1所示)，圖中黑色正三角表示設(shè)計(jì)參考站位置。JLCORS是我國(guó)地區(qū)(省)級(jí)的連續(xù)運(yùn)行的參考站網(wǎng)絡(luò)系統(tǒng)，它將在全省范圍內(nèi)建立的永久性參考站通過(guò)網(wǎng)絡(luò)互聯(lián)，構(gòu)成新一代的網(wǎng)絡(luò)化的大地測(cè)量系統(tǒng)，將GNSS技術(shù)綜合應(yīng)用于吉林省的大地測(cè)量、工程測(cè)量、氣象監(jiān)測(cè)、地震監(jiān)測(cè)、地面沉降監(jiān)測(cè)及城市地理信息系統(tǒng)等領(lǐng)域，2010年12月已經(jīng)建成49個(gè)站并投入運(yùn)行。

2. GNSS后處理軟件

GAMIT/GLOBK是國(guó)際上三大高精度GNSS后處理軟件之一，由美國(guó)麻省理工學(xué)院(MIT)和斯克里普斯海洋研究所(SIO)聯(lián)合開發(fā)，主要應(yīng)用于分析研究地殼形變、高精度GPS測(cè)量數(shù)據(jù)處理等領(lǐng)域。GAMIT軟件處理雙差觀測(cè)量，采用最小二乘算法進(jìn)行參數(shù)估計(jì)。GAMIT軟件主要由以下幾個(gè)模塊構(gòu)成：ARC(軌道積分模塊)、MODEl(組成觀測(cè)方程)、SINCI N(單差自動(dòng)修復(fù)周跳)、DBCI N(雙差自動(dòng)修復(fù)周跳)、CVIEW(人工交互式修復(fù)周跳)、SOI VE(最小二乘解算模塊)、DFMRG(數(shù)據(jù)融合模塊)、FXDRV(生成批處理文件)、GI OBK(運(yùn)用卡爾曼濾波進(jìn)行網(wǎng)平差)等。

圖1　JLCORS基準(zhǔn)站分布

3. 時(shí)間序列噪聲分析方法

時(shí)間序列是動(dòng)態(tài)測(cè)試中被觀測(cè)量在一定的幅值范圍內(nèi)隨時(shí)間出現(xiàn)的一系列隨機(jī)數(shù)據(jù)。頻譜分析和波形分析是最重要和最基本的方法。頻譜分析可求得時(shí)間序列的幅值譜、相位譜、功率譜和各種譜密度等，波形分析也即時(shí)間域分析，與頻譜分析可用傅里葉變換相互轉(zhuǎn)換。功率譜技術(shù)和最大似然法分別在頻率域和時(shí)間域?qū)Υ蟮販y(cè)量手段時(shí)間序列的噪聲分析表明，時(shí)間序列包含的噪聲類型非常復(fù)雜，如白噪聲、冪律噪聲、高斯-馬爾可夫噪聲、帶通濾波噪聲等，可以是僅含有一種噪聲，也可以是多種噪聲的組合。相應(yīng)地用這些噪聲的數(shù)學(xué)模型擬合時(shí)間序列來(lái)求解噪聲參數(shù)。

最大似然估計(jì)(MLE)[9]可以同時(shí)獲得時(shí)間相關(guān)的噪聲模型的所有參數(shù)，經(jīng)過(guò)試驗(yàn)證實(shí)，MLE方法不但可以最大限度地發(fā)現(xiàn)一列隨機(jī)數(shù)據(jù)包含的噪聲，而且能準(zhǔn)確計(jì)算出這些噪聲的功率譜，與譜分析所得的結(jié)果能夠精確符合。本次研究應(yīng)用MLE方法進(jìn)行時(shí)間序列噪聲分析。

(1) 最大似然估計(jì)

最大似然估計(jì)的公式為

(1)

在實(shí)際中，C是一種或多種噪聲統(tǒng)計(jì)模型的組合，可以用下式表示

(2)

式中，σ為噪聲模型的幅度；R為各噪聲模型的協(xié)方差矩陣。

式(2)取自然對(duì)數(shù)后可得

(3)

對(duì)σ求微分可得

(4)

通過(guò)上述變換，可以把多維問題轉(zhuǎn)化為一維問題，在解似然方程中，采用上述單純形法。

(2) 噪聲模型

根據(jù)地球物理觀測(cè)所包含噪聲統(tǒng)計(jì)規(guī)律，構(gòu)造出如下模型，在噪聲分析中擬合時(shí)間序列噪聲。

① 冪律噪聲

因許多地球物理過(guò)程可以用一個(gè)冪律過(guò)程來(lái)描述，GNSS時(shí)間序列噪聲也具有這樣的冪律性質(zhì)，即噪聲功率譜Px(f)與其對(duì)應(yīng)的噪聲頻率f之間存在某種冪次關(guān)系

(5)

(6)

式中，Px為功率譜；f為噪聲頻率；P0、f0為常數(shù)；κ為譜指數(shù)；σpl為噪聲幅度。譜指數(shù)分布范圍通常為(-3，1),穩(wěn)定過(guò)程譜指數(shù)為(-1，1)，不穩(wěn)定過(guò)程譜指數(shù)為(-3，-1)。隨機(jī)過(guò)程分為穩(wěn)定過(guò)程和不穩(wěn)定過(guò)程，不穩(wěn)定過(guò)程在低頻部分通常有較大的功率，而在高頻部分有負(fù)的譜指數(shù)，通常負(fù)譜指數(shù)分布在開區(qū)間(-3，-1)內(nèi)，其中白噪聲的譜指數(shù)κ=0，閃爍噪聲的譜指數(shù)κ=-1、κ=1。類隨機(jī)游走過(guò)程的譜指數(shù)分布在區(qū)間(1，3)內(nèi),如典型的隨機(jī)游走噪聲譜指數(shù)κ=2。

② 高斯-馬爾可夫噪聲

一階高斯-馬爾可夫噪聲統(tǒng)計(jì)模型可用下式表示

(7)

一階高斯馬爾可夫信號(hào)可用下面微分方程表示

(8)

式中，Px為功率譜；f為噪聲頻率；P0、f0為常數(shù)；κ為譜指數(shù)。

③ 帶通噪聲

帶通信號(hào)是周期過(guò)程，在地球物理現(xiàn)象中常常是由季節(jié)變化引起，公式如下

(9)

式中，P為功率譜；f為噪聲頻率；f1為常數(shù)；κ為譜指數(shù)；σbp為噪聲幅度。

二、數(shù)據(jù)處理結(jié)果及噪聲性質(zhì)分析

在獲得JLCORS坐標(biāo)時(shí)間序列時(shí)，使用的數(shù)據(jù)策略是首先在GAMIT部分獲得所有站點(diǎn)的單日松弛解。在計(jì)算中選用IGS站(CHAN、BJFS、KHRJ、AIRA等)和陸態(tài)網(wǎng)站(HLAR、SUIY、JIXN、CHUN等)作為控制點(diǎn)，適當(dāng)緊約束控制站坐標(biāo)而對(duì)JLCORS站給予較松弛的約束。在GLRED模塊，先獲得所有點(diǎn)的時(shí)間序列，對(duì)數(shù)據(jù)結(jié)果進(jìn)行檢查，利用GLOBK將SPOAC給出的全球單日松弛解和計(jì)算所得的區(qū)域單日松弛解進(jìn)行綜合平差計(jì)算，在此基礎(chǔ)上通過(guò)IGS核心站求解相對(duì)于全球參考框架ITRF2008的相似變換參數(shù)，從而獲得ITRF2008下的單日解。最終獲得JLCORS參考站坐標(biāo)時(shí)間序列，分別為E(東向)、N(北向)，U(垂直方向)，如圖2所示。

數(shù)據(jù)處理采用CATS程序[10]，其核心算法是參數(shù)估計(jì)中的最大似然估計(jì)。計(jì)算時(shí)，擬合多個(gè)經(jīng)典噪聲統(tǒng)計(jì)模型到觀測(cè)時(shí)間序列，通過(guò)構(gòu)造和求解似然函數(shù)，獲得最佳的噪聲模型，而且可以同時(shí)求解噪聲模型中的所有待估參數(shù)。計(jì)算是一個(gè)迭代的過(guò)程，為提高處理速度，CATS在計(jì)算過(guò)程中引入1個(gè)標(biāo)量angle，即一個(gè)幫助計(jì)算的輔助量，沒有物理意義，通過(guò)angle的增加計(jì)算出似然函數(shù)的最大似然值，然后計(jì)算噪聲模型中各個(gè)參數(shù)的值。

圖2　JLCORS部分參考站(DGAN、DHUA、EDAO、FJTN)坐標(biāo)時(shí)間序列

通過(guò)對(duì)冪律噪聲、閃爍噪聲、隨機(jī)游走噪聲及它們各自和白噪聲的組合進(jìn)行似然函數(shù)值的計(jì)算，隨機(jī)游走噪聲和白噪聲的組合能取得最大的似然函數(shù)值。表1是觀測(cè)站N-S向坐標(biāo)時(shí)間序列迭代后求得的噪聲組合模型中各個(gè)參數(shù)的估計(jì)值。表中σpl表示利用式(5)的冪律噪聲模型計(jì)算的噪聲幅度，σbp表示利用式(9)的帶通濾波模型計(jì)算的噪聲幅度，SD是各模型計(jì)算結(jié)果的標(biāo)準(zhǔn)差。表中帶通濾波噪聲很多站通過(guò)迭代后，似然函數(shù)不收斂，沒有計(jì)算結(jié)果，用“—”表示，其中有幾個(gè)站有迭代結(jié)果，其標(biāo)準(zhǔn)差也非常大，顯然是錯(cuò)誤的結(jié)果。這樣的情況說(shuō)明序列中不含帶通濾波噪聲，這類噪聲是由于季節(jié)更替產(chǎn)生的，存在于大多地球物理觀測(cè)資料中，在JLCORS參考站中卻沒有解算結(jié)果，所有站坐標(biāo)時(shí)間序列均受白噪聲的影響。由于冪律噪聲模型覆蓋的噪聲類型范圍很廣，所有站結(jié)果可見主要含有冪律噪聲。

表1　JLCORS參考站噪聲參數(shù)分析結(jié)果

用冪律噪聲加白噪聲的組合噪聲模型擬合時(shí)間序列可獲得最大的似然函數(shù)值，同時(shí)獲得每個(gè)參考站所有分向冪律噪聲譜指數(shù)，分布區(qū)間為[-2.32，-1.71]，譜指數(shù)集中在κ=-2附近(如圖3所示)，此時(shí)為冪律噪聲的特例——隨機(jī)游走噪聲。由此可見GNSS時(shí)間序列噪聲占主導(dǎo)地位的應(yīng)該是隨機(jī)游走噪聲。研究地殼運(yùn)動(dòng)的基巖站，隨機(jī)游走噪聲的幅度很小，較短的時(shí)間序列很難監(jiān)測(cè)到，特別是在全球網(wǎng)解算數(shù)據(jù)時(shí)很大的共模誤差會(huì)將其淹沒[11]。本次研究的參考站網(wǎng)為區(qū)域級(jí)別，空間相關(guān)的誤差不會(huì)產(chǎn)生，有利于計(jì)算隨機(jī)游走噪聲的幅度。

圖3　JLCORS東西方向譜指數(shù)和噪聲幅度的關(guān)系

三、結(jié)束語(yǔ)

GNSS觀測(cè)資料結(jié)果解算后在結(jié)果中包含很多誤差信息，如GNSS整個(gè)觀測(cè)系統(tǒng)的誤差、解算過(guò)程的模型誤差、觀測(cè)墩的非構(gòu)造誤差等。在信號(hào)處理過(guò)程中，這些信號(hào)可看作是隨機(jī)噪聲。本文通過(guò)最大似然法進(jìn)行參數(shù)估計(jì)，根據(jù)JLCORS參考站坐標(biāo)時(shí)間序列與經(jīng)典噪聲模型擬合，估計(jì)各個(gè)噪聲模型的參數(shù)，來(lái)確定觀測(cè)結(jié)果的誤差信息可用哪種噪聲模型來(lái)描述。該計(jì)算過(guò)程使用噪聲分析軟件CATS來(lái)實(shí)現(xiàn)，經(jīng)過(guò)分析表明，所有參考站均可用冪律噪聲和白噪聲的組合模型來(lái)描述，根據(jù)譜指數(shù)分布區(qū)間，冪律噪聲中的特殊情況之一的隨機(jī)游走噪聲占主導(dǎo)地位。也就是說(shuō)JLCORS參考站坐標(biāo)時(shí)間序列的誤差信息可用隨機(jī)游走噪聲和白噪聲的組合來(lái)描述。

[1]王敏, 沈正康, 董大南. 非構(gòu)造形變對(duì) GPS 連續(xù)站位置時(shí)間序列的影響和修正[J]. 地球物理學(xué)報(bào), 2005, 48(5): 1045-1052.

[2]田云鋒, 沈正康. GPS 坐標(biāo)時(shí)間序列中非構(gòu)造噪聲的剔除方法研究進(jìn)展[J]. 地震學(xué)報(bào), 2009, 31(1): 68-81.

[3]黃立人, 符養(yǎng). GPS 連續(xù)觀測(cè)站的噪聲分析[J]. 地震學(xué)報(bào), 2007, 29(2): 197-202.

[4]LANGBEIN J. Noise in Two-color Electronic Distance Meter Measurements Revisited[J]. Journal of Geophysical Research: Solid Earth (1978—2012), 2004, 109(B4).DOI:10.1029/2003JB002819.

[5]LANGBEIN J, JOHNSON H. Correlated Errors in Geodetic Time Series: Implications for Time-dependent Deformation[J]. Journal of Geophysical Research: Solid Earth (1978—2012), 1997, 102(B1): 591-603.

[6]MAO A, HARRISON C G, DIXON T H. Noise in GPS Coordinate Time Series[J]. Journal of Geophysical Research: Solid Earth (1978—2012), 1999, 104(B2): 797-816.

[7]杜向鋒, 張興福, 張永毅. CORS 測(cè)量成果轉(zhuǎn)換的一步法及其精度分析[J]. 測(cè)繪通報(bào), 2015(7): 23-26.

[8]鄭立平. 港口鎮(zhèn)單參考站 CORS 系統(tǒng)建設(shè)及應(yīng)用[J]. 測(cè)繪通報(bào), 2011(1): 43-45.

[9]莫惠棟. 最大似然法及其應(yīng)用[J]. 遺傳, 1984, 6(5): 42-48.

[10]WILLIAMS S D. CATS: GPS Coordinate Time Series Analysis Software[J]. GPS solutions, 2008, 12(2): 147-153.

[11]BOS M, FERNANDES R, WILLIAMS S, et al. Fast Error Analysis of Continuous GPS Observations[J]. Journal of Geodesy, 2008, 82(3): 157-166.

Analysis of the Noise of Coordinate Time Series from JLCORS Network Reference Stations

LIU Junqing,XIAO Huifeng,ZHANG Chenxia,DING Guang

劉俊清，肖輝峰，張晨俠,等.JLCORS網(wǎng)絡(luò)參考站坐標(biāo)時(shí)間序列噪聲類型分析[J].測(cè)繪通報(bào)，2016(10)：1-5.DOI:10.13474/j.cnki.11-2246.2016.0316.

2015-12-15

地震科技星火計(jì)劃(XH14016Y)

劉俊清(1977—)，男，博士生，高級(jí)工程師，主要從事地震火山監(jiān)測(cè)預(yù)報(bào)研究。E-mail: lunwen_98@126.com

P228

B

0494-0911(2016)10-0001-05