楊 博,占 偉,高艷龍,韓月萍
(中國(guó)地震局第一監(jiān)測(cè)中心,天津 300180)
GNSS垂向形變場(chǎng)信息提取的嘗試
楊 博,占 偉,高艷龍,韓月萍
(中國(guó)地震局第一監(jiān)測(cè)中心,天津 300180)
一定條件下GNSS流動(dòng)復(fù)測(cè)資料可為我們揭示地形變提供必要的依托,但由于垂向分量觀測(cè)誤差較大而沒有被恰當(dāng)?shù)乩?。利用華北地區(qū)1999—2007年的4期資料,從形變場(chǎng)的連續(xù)變化角度出發(fā),提出了利用多核函數(shù)進(jìn)行解析與濾波的方法,理論與結(jié)果表明:①可獲得垂直運(yùn)動(dòng)場(chǎng)的數(shù)學(xué)解析;②可實(shí)現(xiàn)對(duì)垂向運(yùn)動(dòng)場(chǎng)不同層次的濾波與信息分離;③當(dāng)流動(dòng)測(cè)站具有一定的空間密度、一定數(shù)量的復(fù)測(cè)結(jié)果和時(shí)間長(zhǎng)度則可獲得較可靠的形變信息。最后利用華北地區(qū)GNSS流動(dòng)觀測(cè)資料給出了實(shí)際算例與初步分析結(jié)果。
GNSS;垂向運(yùn)動(dòng);多核函數(shù);濾波與解析
在地球動(dòng)力學(xué)研究領(lǐng)域里,以GNSS觀測(cè)資料為依托的研究已是其重要的一個(gè)分支;具體到地殼形變研究而言,利用GNSS資料人們有望從空間與時(shí)間的角度觀察地殼三維的運(yùn)動(dòng)及形變的詳細(xì)過程。因此,多角度的研究也相應(yīng)地展開[1-12],但目前基本都是處在二維水平形變的研究層面上。這主要是因?yàn)榱鲃?dòng)GNSS垂向觀測(cè)誤差較大、影響的因素也較多,所以垂向觀測(cè)分量在形變分析中一直沒有得到應(yīng)用。然而,隨著復(fù)測(cè)資料的積累,垂向運(yùn)動(dòng)處理的信噪比也在不斷的提高,這是一個(gè)不爭(zhēng)的事實(shí)。為了更加較有效地利用垂向形變信息,本文試圖從運(yùn)動(dòng)“觀測(cè)值”入手,在形變場(chǎng)連續(xù)變化的前提下,通過一定的數(shù)理方法,實(shí)現(xiàn)垂直形變場(chǎng)的空間解析與濾波,以更快更好地利用垂向運(yùn)動(dòng)信息并服務(wù)于地殼形變研究、地表監(jiān)測(cè)等領(lǐng)域。另一方面,大規(guī)模地開展精密水準(zhǔn)測(cè)量在較短的時(shí)段內(nèi)已不太可能,就華北而言,已近20年沒有進(jìn)行過全面的復(fù)測(cè)。所以,開發(fā)利用其他手段的監(jiān)測(cè)已是非常必要了,GNSS或許是當(dāng)今最切實(shí)際的一種選擇。
多核函數(shù)法是Hardy在1971年提出的一種數(shù)
式中sj(x,y,xj,yj)為核函數(shù),(xj,yj)為核點(diǎn)的坐標(biāo),cj為待定系數(shù)。
一般說核函數(shù)是一個(gè)非常簡(jiǎn)單的曲面形態(tài),通常為“缽”形或“鐘”形曲面函數(shù)。盡管哪一種形態(tài)的核函數(shù)更優(yōu)尚無(wú)非常明確的定論,但在此后的應(yīng)用研究中,該法已被廣泛地用于地殼形變研究[13-14]、DEM 內(nèi)插[15]和地球重力場(chǎng)元內(nèi)插[16],并被證明是非常成功和有效的。
然而,多核函數(shù)在具體應(yīng)用中的效果如何仍存在選擇問題。對(duì)于不同的描述對(duì)象可能會(huì)有不同的選擇,即使對(duì)相同的對(duì)象也并非是唯一選擇,甚至很難區(qū)分它們之間的優(yōu)劣。所以,在實(shí)際應(yīng)用中主要取決于使用者對(duì)多核函數(shù)及描述對(duì)象的理解。一般取決于3個(gè)方面:①核函數(shù)的選擇;②核函數(shù)中的參數(shù)選擇;③核函數(shù)的核點(diǎn)選擇。由于我們所描述的對(duì)象是較大尺度的地殼形變場(chǎng),而地殼介質(zhì)由于其各向非同性、活動(dòng)構(gòu)造帶的存在及應(yīng)力大小與方向的空間變化等,使得形變場(chǎng)具有廣譜性。這就是說,我們所描述的對(duì)象是很復(fù)雜的。據(jù)此,這里所選用的核函數(shù)為“缽”形函數(shù),具體為:值逼近方法,其基本思想是,任何數(shù)學(xué)表面和任何不規(guī)則的圓滑表面,總可以通過一系列有規(guī)則的數(shù)學(xué)表面的合成,以任意精度逼近。其數(shù)學(xué)表達(dá)式為:
式中:dj為球面上兩點(diǎn)間的大地線長(zhǎng)度(以km為單位),ST= (s1,s2,…,snx),CT= (c1,c2,…,cnx)。
在GNSS計(jì)算中通常是以ITRF(或WGS84)作為參考框架獲取速度場(chǎng),因此測(cè)站j的位置(λj,φj)、運(yùn)動(dòng)速率νu(λj,φj)及其相應(yīng)的誤差mu(λj,φj)都是已知的。在用多核函數(shù)進(jìn)行數(shù)值逼近時(shí)應(yīng)將所有測(cè)站位置作為核函數(shù)的核點(diǎn)位置,因此可列出n個(gè)方程(測(cè)站的個(gè)數(shù)),故可求解n個(gè)未知數(shù)(C)。這時(shí)即可獲得垂向運(yùn)動(dòng)的解析式fu(λ,φ),同理也可獲得誤差解析式fm(λ,φ)。即
由于上述數(shù)值逼近所得到的是垂向運(yùn)動(dòng)觀測(cè)值的解析結(jié)果,它不僅包含運(yùn)動(dòng)信息同時(shí)也包含誤差干擾(運(yùn)動(dòng)觀測(cè)值中所帶來(lái)的),所以由此得到的運(yùn)動(dòng)結(jié)果實(shí)際上是各種成分的綜合結(jié)果。因此,進(jìn)行濾波與信息分離是必要的。多核函數(shù)不但可以進(jìn)行數(shù)值逼近,而且也具有濾波的功能,故這里進(jìn)行空間濾波同樣利用多核函數(shù)法。進(jìn)行濾波計(jì)算時(shí),濾波效果或?yàn)V波層次則取決于核點(diǎn)的選用與分布。從性質(zhì)上看,核點(diǎn)的間隔越?。ú坏眯∮跍y(cè)站間的平均間隔),由此所獲得的運(yùn)動(dòng)信息所包含的頻譜帶寬越寬,可含蓋從低頻至高頻的范圍,反之帶寬越窄,也越趨于低頻?,F(xiàn)假定垂向運(yùn)動(dòng)分量的濾波函數(shù)為:
式中:AT=(a1,…,anx)為待定系數(shù)矩陣;ST= (s1,…,snx)為 核 函 數(shù) 矩 陣,si(λ,φ,λi,φi)= exp(-(di/780)2.8)(di為球面上兩點(diǎn)間的大地線長(zhǎng)度,以km 為單位),(λi,φi)為核點(diǎn)位置坐標(biāo)。
在實(shí)際計(jì)算時(shí),依據(jù)解析式(4)首先以相同的網(wǎng)格化形式(等步長(zhǎng))獲取垂向分量和相應(yīng)誤差的數(shù)據(jù),并以垂向分量作為“觀測(cè)值”。然后,根據(jù)最小二乘法求解(5)式中的待定系數(shù),即通過誤差方程:
其中N=DTPL,W=DTPD,L為相應(yīng)分量觀測(cè)值的矩陣,P 為相應(yīng)分量觀測(cè)值的權(quán)矩陣(pi=1/m2i)。當(dāng)X求得后(即(5)式中的A),(5)式則變?yōu)榻?jīng)濾波后垂向運(yùn)動(dòng)計(jì)算的函數(shù)解析式,據(jù)此可獲得任意位置垂向運(yùn)動(dòng)濾波值。
根據(jù)誤差傳播律,可進(jìn)行濾波值的誤差估計(jì)。由(6)式等可計(jì)算單位權(quán)方差估值:
式中nx為核點(diǎn)的個(gè)數(shù),n為測(cè)站數(shù)量,且大于nx(否則無(wú)法給出精度評(píng)定)。垂向運(yùn)動(dòng)方差則為:
考慮到信息的全面性,還可進(jìn)一步給出梯度在空間上的變化結(jié)果。由于以上已經(jīng)獲得了運(yùn)動(dòng)場(chǎng)的解析式,所以梯度計(jì)算則變得非常方便了。根據(jù)梯度的定義,有:
據(jù)此就可很容易地給出梯度及其在空間的變化結(jié)果,并服務(wù)于垂直形變的全面分析。
華北上地殼最突出的現(xiàn)今動(dòng)力學(xué)表現(xiàn)是新生代盆地的發(fā)育,而且盆地基底面的起伏與莫霍面的起伏呈上凹下凸型的鏡像關(guān)系。這實(shí)際上表明深部活動(dòng)對(duì)淺部構(gòu)造的制約活動(dòng),其垂直活動(dòng)特征為山區(qū)隆升,平原下降,并稱之為常態(tài)活動(dòng)。華北地形地貌的空間特征如圖1所示:由北京-石家莊-鄭州連線以西以北地區(qū)為山區(qū)或海拔較高的地區(qū);其次是山東地區(qū),其中較大一部分山區(qū)或丘陵地區(qū),在其南部的安徽地區(qū)也有一些具有一定海拔高度的丘陵地貌;除此之外,基本上為海拔較低的平原地區(qū)(京津及以南的河北、河南和江蘇等地區(qū))。由于華北地區(qū)是我國(guó)地震活動(dòng)較為強(qiáng)烈的地區(qū)之一,因此中國(guó)地殼運(yùn)動(dòng)網(wǎng)絡(luò)在此布設(shè)了較密集的GNSS流動(dòng)測(cè)站(測(cè)站分布為圖2中的黑點(diǎn)所示),以用于監(jiān)測(cè)該地區(qū)地殼形變及其動(dòng)態(tài)變化。
圖1 華北地形地貌空間特征略圖
GNSS前期數(shù)據(jù)處理主要用較新版本GAMIT/GLOBK/QOCA軟件完成,其步驟可概括為:①利用GAMIT軟件獲得測(cè)站坐標(biāo)和衛(wèi)星軌道的單日松弛解;②利用GLOBK軟件將獲得的每天單日松弛解和IGS數(shù)據(jù)中心SOPAC產(chǎn)出的全球IGS跟蹤站有關(guān)的單日松弛解合并;③利用QOCA軟件綜合所有的單日松弛解計(jì)算各測(cè)站的三維運(yùn)動(dòng)速度。在利用QOCA軟件求解測(cè)站速度時(shí),采用ITRF05地球框架的速度解作為約束,聯(lián)合解算了1999—2007年(1999、2001、2004和2007年4期觀測(cè)資料)這一時(shí)段各測(cè)站包括垂向的三維運(yùn)動(dòng)速度。在由GNSS所獲得的垂向運(yùn)動(dòng)速度的基礎(chǔ)上,利用文中所述的方法進(jìn)行進(jìn)一步處理。其基本步驟是,首先依據(jù)(3)、(4)式對(duì)研究區(qū)相對(duì)運(yùn)動(dòng)速度及誤差結(jié)果進(jìn)行數(shù)據(jù)逼近,以獲得解析式;然后以此基礎(chǔ)以經(jīng)、緯步長(zhǎng)分別為35km間隔計(jì)算網(wǎng)格點(diǎn)上的垂向運(yùn)動(dòng)結(jié)果;接下來(lái)再利用(5)~(8)式進(jìn)行空間濾波計(jì)算。濾波計(jì)算時(shí)選用的核點(diǎn)其經(jīng)、緯間隔均為150km,結(jié)果見圖2。
圖2所顯示的是華北地區(qū)1999—2007年垂向運(yùn)動(dòng)的趨勢(shì)性結(jié)果,由這樣的結(jié)果我們可獲得第一印象是華北地區(qū)近十來(lái)年的運(yùn)動(dòng)形態(tài)在性質(zhì)上與地形地貌形態(tài)是一致的,零值線的展布基本上體現(xiàn)了山區(qū)與平原的分界線;它展現(xiàn)了山區(qū)、丘陵等海拔較高的地區(qū)為相對(duì)隆升區(qū)(圖2中的藍(lán)色區(qū)域),海拔較低的平原地區(qū)為相對(duì)下沉區(qū)(圖2中的紅色區(qū)域)的特征。換句話說,華北地區(qū)當(dāng)前的垂向運(yùn)動(dòng)方式仍處于繼承性運(yùn)動(dòng)方式。但就空間的細(xì)部變化來(lái)看,不同的地區(qū)仍有較大的差別。華北地區(qū)西部的鄂爾多斯及周邊地區(qū)整體雖表現(xiàn)為上升,平均上升接近3mm/a,但內(nèi)部的差異變化并不大,變化范圍0~5mm/a;此外,由于等值線分布寬緩,這表明其梯度變化是很小的,體現(xiàn)了內(nèi)部運(yùn)動(dòng)具有很好的一致性;從地震的孕育角度看,內(nèi)部差異運(yùn)動(dòng)較小的地域或塊體難以儲(chǔ)存較高的能量,因此也難以發(fā)生較大的地震,該區(qū)的地震歷史也表明了這一點(diǎn),同時(shí)也說明這一狀態(tài)還將持續(xù)下去。鄂爾多斯東邊緣的山西地域是華北地區(qū)隆升最突出的地域,除總體表現(xiàn)為上升外,最大上升量也位于其中,達(dá)10mm/a,其地點(diǎn)在太原以東地區(qū)。在山西境內(nèi)隆升垂向運(yùn)動(dòng)變化范圍為0~10mm/a,這表明山西地域的垂直形變是比較突出的,但并非均勻,主要集中在太原的周圍地區(qū)。由于最小隆升位于太原以西地區(qū),因此與其東側(cè)最大隆升相比形成了鮮明的對(duì)照,說明該較小區(qū)域是能量相對(duì)積累最快的區(qū)域。除此之外,山西地域垂向運(yùn)動(dòng)的另一個(gè)特征是山西斷陷帶相對(duì)較小,其東部相對(duì)較大,具有區(qū)域上的繼承性質(zhì)。華北北端燕山帶的整體隆升表現(xiàn)為“西小東大”,變化范圍為2~6mm/a,形變相對(duì)均勻。山東的山區(qū)和丘陵地域基本上為隆升運(yùn)動(dòng),但變化不均勻,大體以郯廬帶為界,兩側(cè)上升,具有一定的對(duì)稱形態(tài),西側(cè)最大上升為6mm/a,東側(cè)為8mm/a,并存在著變形。其南的安徽地區(qū)也以隆升運(yùn)動(dòng)為主,最大達(dá)6mm/a,但等值線分布較為寬緩。以上所述的隆升運(yùn)動(dòng),基本上體現(xiàn)的是近十年來(lái)的構(gòu)造運(yùn)動(dòng)。然而,華北平原、江蘇等地區(qū)運(yùn)動(dòng)發(fā)生了很大的變化,它們均處在下沉的狀態(tài)之中。就性質(zhì)而言,更多的屬于地面下沉,尤其是石家莊以東、北京以南、濟(jì)南以西等所圍地區(qū),明顯表現(xiàn)為一個(gè)巨型的沉降“漏斗”,其中心沉降接近50mm/a。就其形成機(jī)理而言,基本上是過量開采地下水所致,并非屬于構(gòu)造運(yùn)動(dòng)的反映。河南和江蘇地區(qū)等下沉量較小,一般在5mm/a以內(nèi),形變相對(duì)也較小。由于其“漏斗”的特征不明顯,所以構(gòu)造運(yùn)動(dòng)的比重相對(duì)較大些。
圖2 華北地區(qū)1999—2007年垂向運(yùn)動(dòng)的趨勢(shì)性結(jié)果(單位:mm/a)
如果不考慮沉降漏斗這一空間范圍,圖2所示結(jié)果基本上表征了華北近十年來(lái)具有構(gòu)造活動(dòng)含義的垂向形變場(chǎng)。因此,它告訴我們當(dāng)前華北地區(qū)的垂向運(yùn)動(dòng)是繼承性的。但在繼承性運(yùn)動(dòng)中,或許存在著變形偏大的運(yùn)動(dòng)空間,如太原周圍地區(qū),故應(yīng)給予關(guān)注。事實(shí)上,華北地區(qū)十幾年來(lái)并未發(fā)生強(qiáng)烈地震;根據(jù)歷史地震,一方面說明華北存在缺震現(xiàn)象;另一方面也可能說明缺震是有一定背景的。該區(qū)的繼承性運(yùn)動(dòng)表明能量還處在積累之中,但對(duì)形變較大的部位,及無(wú)形變且歷史上又是多震的部位不能掉以輕心。
綜上所述,在一定條件下GNSS的垂向運(yùn)動(dòng)信息是可以利用的,如果在GNSS的前期與后期數(shù)據(jù)處理進(jìn)一步開展工作,可以獲得服務(wù)于不同領(lǐng)域而又有效的形變信息。本文主要圍繞后期數(shù)據(jù)方法而進(jìn)行一些探討和嘗試。通過多核函數(shù)法實(shí)現(xiàn)了運(yùn)動(dòng)場(chǎng)的解析和濾波計(jì)算,同時(shí)也給出了運(yùn)動(dòng)場(chǎng)空間濾波結(jié)果的嚴(yán)密精度評(píng)定算式。實(shí)際結(jié)果表明,由此所獲得的形變信息是有參考價(jià)值的。這在當(dāng)前利用精密水準(zhǔn)進(jìn)行較大空間尺度垂直形變場(chǎng)研究難以維系的情況下,可以說是希望。然而,本方法更加恰當(dāng)?shù)氖褂眠€需結(jié)合實(shí)際情況、研究目標(biāo)、及測(cè)站的分布等。核函數(shù)點(diǎn)的選擇與濾波和信息分離有密切的關(guān)系,一般地說,核函數(shù)點(diǎn)間隔越大則信息越趨于低頻,反之亦然。
通過對(duì)華北地區(qū)GNSS復(fù)測(cè)資料的處理與垂向形變信息的提取,使我們獲得的最基本看法是,華北地區(qū)的構(gòu)造運(yùn)動(dòng)在1999—2007年期間以繼承性運(yùn)動(dòng)為主,即山區(qū)或海拔較高的地區(qū)為隆升,隆升量級(jí)絕大部分在4mm/a以內(nèi),極個(gè)別區(qū)域達(dá)8~10mm/a,平原地區(qū)為下降,可辨別的構(gòu)造運(yùn)動(dòng)約在5mm/a以內(nèi),華北平原的下降主要屬于過量開采地下流體造成的,其下沉量最大達(dá)50mm/a。
致謝:本項(xiàng)研究得到了楊國(guó)華研究員的具體指導(dǎo),在此表示感謝。
[1] 楊國(guó)華,張風(fēng)霜,武艷強(qiáng),等.利用 GPS連續(xù)觀測(cè)資料進(jìn)行強(qiáng)震危險(xiǎn)性預(yù)測(cè)的探討[J].地震,2008,28(1):33-39.
[2] 黃立人,王敏.中國(guó)大陸構(gòu)造塊體的現(xiàn)今活動(dòng)和變形[J].地震地質(zhì),2003,25(1):23-32.
[3] 楊國(guó)華,李延興,韓月萍,等.由 GPS觀測(cè)結(jié)果推導(dǎo)中國(guó)大陸現(xiàn)今水平應(yīng)變場(chǎng)[J].地震學(xué)報(bào),2002,24(4):337-347.
[4] 江在森,馬宗晉,張希.GPS初步結(jié)果揭示的中國(guó)大陸水平應(yīng)變場(chǎng)與構(gòu)造變形[J].地球物理學(xué)報(bào),2003,46(3):352-358.
[5] 楊國(guó)華,江在森,王敏.印尼地震對(duì)我國(guó)川滇地區(qū)地殼水平活動(dòng)的影響[J].大地測(cè)量與地球動(dòng)力學(xué),2006,26(1):9-14.
[6] 楊博,周偉,陳阜超,等.GPS連續(xù)站水平位置時(shí)間序列共模白噪聲識(shí)別與估計(jì)的歐拉-濾波法[J].山東科技大學(xué)學(xué)報(bào)(自然科學(xué)版),2010,29(3):26-31.
[7] 楊國(guó)華,楊博,張風(fēng)霜,等.汶川地震對(duì)華北地區(qū)水平形變場(chǎng)影響及有關(guān)含義的討論[J].地震,2009,29(1):77-83.
[8] 楊國(guó)華,韓月萍,楊博,等.華北地區(qū)最近幾年水平形變場(chǎng)的含義與討論[J].山東科技大學(xué)學(xué)報(bào)(自然科學(xué)版),2009,28(6):1-8.
[9] 王敏,沈正康,牛之俊,等.現(xiàn)今中國(guó)大陸地殼運(yùn)動(dòng)與活動(dòng)塊體模型[J].中國(guó)科學(xué)(D輯),2003,33(增):21-32.
[10] 王琪,丁國(guó)瑜,喬學(xué)軍,等.天山現(xiàn)今地殼快速縮短與南北地塊的相對(duì)運(yùn)動(dòng)[J].科學(xué)通報(bào),2000,45(14):1543-1547.
[11] 楊博,陳阜超,周偉,等.GPS連續(xù)站水平位置坐標(biāo)時(shí)間序列插值的一種新方法[J].華北地震科學(xué),2009,28(2):22-27.
[12] 楊國(guó)華,江在森,武艷強(qiáng),等.中國(guó)大陸整體無(wú)凈旋轉(zhuǎn)基準(zhǔn)及其應(yīng)用[J].大地測(cè)量與地球動(dòng)力學(xué),2005,29(4):6-10.
[13] 黃立人,劉天奎.用速度面擬合法研究華北部分地區(qū)的現(xiàn)今地殼垂直運(yùn)動(dòng)[J].地殼形變與地震,1989,9(3):6-12.
[14] 楊國(guó)華,黃立人.速度面擬合法中多面函數(shù)幾個(gè)特性的初步數(shù)值研究[J].地殼形變與地震,1990,10(4):70-82.
[15] 周光文,黃莜容.一種改進(jìn)的多面函數(shù)法[J].測(cè)繪通報(bào),1996,(6):6-8.
[16] 陶本藻.GPS水準(zhǔn)似大地水準(zhǔn)面擬合和正常高計(jì)算[J].測(cè)繪通報(bào),1992,(4):14-18.
[17] 陳阜超,紀(jì)靜,塔拉,等.京津水準(zhǔn)復(fù)測(cè)與垂直形變特征[J].華北地震科學(xué),2011,(2):31-34.
Some Attempt about Distilling Information of Vertical Deformation Field from GNSS Data
YANG Bo,ZHAN Wei,GAO Yan-long,HAN Yue-ping
(First Crust Monitoring and Application Center,CEA,Tianjin 300180,China)
In certain condition,GNSS mobile observations can provide essential support to reveal crustal deformation,but they were not used properly because of large observation error in vertical component.From the perspective of continuous variation of deformation field,a method by using multi-kernel function to filtering and analysis is proposed with 4observations(1999—2007)in the GPS network of North China,theory and results show as follows.Firstly,mathematical analysis of vertical motion field can be obtained.Secondly,different levels of filtering and information separation of vertical deformation field can be realized.Thirdly,if the space density of mobile stations,quantity and time length of repetition measurement can reach a certain level,reliable deformation information can be achieved.Finally,apractical example and preliminary analysis results with GNSS mobile observations in North China were given.
GNSS;vertical movement;multi-kernel function;filtering and analysis
P315.72
A
1003-1375(2012)01-0001-05
2011-07-15
地震行業(yè)科研專項(xiàng)重大項(xiàng)目(200908029);地震行業(yè)科研專項(xiàng)(201008007)
楊博(1985-),男(漢族),天津人,助理工程師,主要從事GNSS數(shù)據(jù)處理與分析工作.E-mail:boyyang1234@qq.com