薄萬舉 鄭智江 武艷強 張立成 張雙喜 王同慶 暢 柳

1 中國地震局第一監(jiān)測中心，天津市耐火路7號，300180

精密水準(zhǔn)是用于測量地表點位高程的傳統(tǒng)大地測量手段之一[1-2]；GNSS測量則是近30 a來迅速發(fā)展起來的新技術(shù)，能快速給出地面點的三維坐標(biāo)。研究表明，對于同一個點，二者給出的垂向坐標(biāo)存在差異，但屬于近似不變的系統(tǒng)性差異，其對點位垂向變化量計算的影響可以忽略。因此，通過分析不同時段觀測結(jié)果的變化，二者都能較好地給出地表點的垂直形變量[3]。但實踐中發(fā)現(xiàn)，不同來源的資料、不同技術(shù)手段得到的資料、相同資料用不同處理方法給出的結(jié)果往往不同，有時差異很大，甚至完全相反。對比以精密水準(zhǔn)給出的全國垂直形變矢量圖和GNSS給出的垂直形變矢量圖發(fā)現(xiàn)，兩者在天山構(gòu)造帶及塔里木盆地區(qū)域差異最大，在天山構(gòu)造帶甚至完全相反。本文針對這一差異進(jìn)行分析與研究，探討其中可能存在的問題，分析形成差異的原因。

1 以精密水準(zhǔn)為主給出的垂直形變矢量圖

圖1是以精密水準(zhǔn)給出的天山構(gòu)造帶及塔里木盆地區(qū)域垂直形變矢量圖(割圓錐投影)。所用資料為2期全國精密水準(zhǔn)資料，第1期中心時刻為1980年，第2期中心時刻為2010年，獲得該期間全國垂直形變圖后，截取天山構(gòu)造帶區(qū)域部分，即為圖1。為了抑制精密水準(zhǔn)隨路線長度增加而導(dǎo)致的誤差累積的影響，優(yōu)選全國57個GNSS連續(xù)站垂向運動速率作為約束(均為與水準(zhǔn)點有2期以上聯(lián)測成果的GNSS連續(xù)站，其分布見圖2)，利用赫爾默特平差方法[4]得到各個水準(zhǔn)點的垂向運動速率。全部GNSS站點坐標(biāo)運動速率是在全球GNSS參考框架下統(tǒng)一解算的，用57個GNSS連續(xù)站的垂向速率對全國水準(zhǔn)網(wǎng)的速率平差進(jìn)行約束，因GNSS站的垂向速率統(tǒng)計精度相對均勻且高于考慮誤差傳播和累積后由單一水準(zhǔn)測量得到的垂向運動速率的精度，約束的結(jié)果相當(dāng)于統(tǒng)一了精密水準(zhǔn)與GNSS的垂向速率基準(zhǔn)。當(dāng)然這只有在地殼垂向速率為線性的前提下才是嚴(yán)密的，實際上地殼垂直運動速率的變化是非線性的，水準(zhǔn)與GNSS資料覆蓋的時段又不完全吻合，因此二者速率基準(zhǔn)的統(tǒng)一是相對的和近似的，尚存在一定的誤差。

圖1 以精密水準(zhǔn)為主給出的天山構(gòu)造帶及 塔里木盆地區(qū)域垂直形變速率矢量圖 (1980～2010年)Fig.1 Vector diagram of vertical deformation rate in Tianshan tectonic belt and Tarim basin mainly based on precision leveling (1980-2010)

圖2 精密水準(zhǔn)路線及2期57個聯(lián)測 GNSS連續(xù)站分布示意圖Fig.2 Schematic diagram of precision leveling routes and 57 GNSS stations with 2 times of joint measurement

以精密水準(zhǔn)為主給出的垂直形變圖顯示，該區(qū)域的水準(zhǔn)點總體上以相對沉降運動為主。

由圖2可見，西部符合條件的GNSS聯(lián)測點較少，本文研究范圍內(nèi)(天山構(gòu)造帶附近)只有2個點，擴大到塔里木盆地有4個點，這4個聯(lián)測點的情況見表1。

表1 天山構(gòu)造帶及塔里木盆地區(qū)域4個GNSS 聯(lián)測點的相關(guān)參數(shù)

對比表1和圖1可見，XJHT、XJSS、XJYC站的垂向速率都對附近的水準(zhǔn)速率起到了很好的約束作用；而XJKL站的GNSS垂向速率為7.0 mm/a，在圖1中似乎對附近的水準(zhǔn)速率失去了約束作用，經(jīng)核查發(fā)現(xiàn)，該站的垂向運動速率確實高。7.0 mm/a的運動速率是真實的，水準(zhǔn)聯(lián)測數(shù)據(jù)也恰好證明了這一點。大面積形變重點研究空間上相對低頻的垂直形變信息，在GNSS垂向形變圖中通過統(tǒng)一的判別標(biāo)準(zhǔn)將XJKL點作為奇異點進(jìn)行剔除，因此在后面的GNSS形變圖中未體現(xiàn)該點突出的垂直運動。但該點恰好有2期與水準(zhǔn)聯(lián)測的成果，資料完整，考慮該區(qū)域缺少符合條件的聯(lián)測點，故在聯(lián)合平差中仍將其作為聯(lián)測點使用，同樣起到了很好地約束于統(tǒng)一速率基準(zhǔn)的作用。

2 GNSS給出的垂直形變速率矢量圖

GNSS垂向分量精度較低，但與精密水準(zhǔn)相比，不存在長距離誤差累積和傳遞的問題，各GNSS點之間具有一定的相對獨立性，當(dāng)GNSS資料足夠豐富時，給出的垂向形變信息具有一定的可信性。中國大陸近20 a來積累了比較豐富的GNSS資料，筆者嘗試用GNSS資料給出一張中國大陸GNSS垂直形變速率矢量圖，截取其中天山構(gòu)造帶及塔里木盆地區(qū)域得到圖3。

圖3 GNSS給出的天山構(gòu)造帶及塔里木盆地區(qū)域 垂直形變速率矢量圖(1999～2018年)Fig.3 Vector diagram of vertical deformation rate in Tianshan tectonic belt and Tarim basin based on GNSS (1999-2018)

3 差異性分析

對比圖1和圖3可以看出，精密水準(zhǔn)給出的垂向形變幾乎是整體沉降，GNSS給出的垂向形變幾乎是整體抬升，即使所跨時段稍有不同，也不應(yīng)該有如此明顯的差異。

首先考慮兩者測點空間分布差異的影響。GNSS測點要求盡量選在穩(wěn)定的地方，且希望能同時收到盡可能多的衛(wèi)星信號，一般要遠(yuǎn)離斷層、河流、低洼等不穩(wěn)定和衛(wèi)星高度角不滿足要求的地方，因此，GNSS觀測站點多選在山頭等視線開闊、海拔相對較高的地方；而水準(zhǔn)路線的布設(shè)則要重點考慮便于觀測的問題，往往沿道路兩邊選點，尤其在山區(qū)，遇到高山能繞開必須繞開，往往要從高山之間的峽谷穿過。這就使GNSS和水準(zhǔn)的重合點非常少，僅有少量幾個經(jīng)過聯(lián)測并用于平差約束的GNSS連續(xù)站可作為二者的重合點。但就是這少量重合點，也可將二者的垂直形變速率基準(zhǔn)統(tǒng)一起來，不應(yīng)該有如此大的差異。

將二者的矢量合在一起，直接繪到同一張圖上(圖4，墨卡托投影,單位mm/a)則可以看出，二者雖然總體上升降相互矛盾，但具體的點位分布卻相互錯開。正如前文所述，水準(zhǔn)點多為沉降，處于相對低的位置；GNSS點多為隆升，處于相對高的位置。高低之差形成的原因主要是垂向構(gòu)造活動存在差異所致，因此正常情況下點位垂直運動速率存在差異也在情理之中。仔細(xì)分析可見，圖4中仍有少量相互矛盾的地方，考慮到觀測誤差、圖1和圖3觀測數(shù)據(jù)所跨時段的差異以及地殼構(gòu)造活動在空間的非均勻性和時間上的非線性等多因素影響，存在少量相互矛盾的地方是可以理解的。

圖4 天山構(gòu)造區(qū)垂直形變速率矢量圖 (圖1和圖3資料合并繪制，墨氏投影)Fig.4 Vector diagram of vertical deformation rate in Tianshan tectonic region(combined drawing of data in Fig.1 and Fig.3,Mercator projection)

利用圖4數(shù)據(jù)，采用克里金插值算法繪制等值線圖(見圖5，單位mm/a)。圖5顯示，虛線勾畫出的沉降區(qū)恰好為海拔相對低的區(qū)域，實線勾畫出的隆起區(qū)多為山區(qū)，與構(gòu)造地貌吻合較好，在一定程度上較好地反映了該區(qū)域在該時段垂向構(gòu)造活動的空間差異。

圖5 天山構(gòu)造區(qū)垂直形變速率等值線圖 (圖1和圖3資料合并繪制，墨氏投影)Fig.5 Contour map of vertical deformation rate in Tianshan tectonic region (combined drawing of data in Fig.1 and Fig.3,Mercator projection)

4 初步融合結(jié)果及解釋

圖5給出的結(jié)果是用純數(shù)學(xué)方法通過空間內(nèi)插得到的等值線。目前常用的通過內(nèi)插數(shù)據(jù)繪制等值線的方法有多種，如最小二乘配置方法[2,5]、面函數(shù)法[6]等。不同方法給出的等值線分布具有很大差異。分析認(rèn)為，不同的方法和模型各有利弊，所用的數(shù)學(xué)模型不同，實際上對應(yīng)著不同的物理意義，有的偏向于低頻形變信息，會濾掉很多高頻形變成分；有的偏向中高頻信息，則掩蓋了低頻形變信息的變化。因此需要根據(jù)實際情況選擇相對匹配的方法或方法組合，以得到更為理想的結(jié)果。對于圖4分布的垂向矢量，有的兩點間距離很近，但高差相差很大，垂向速率相反；有的兩點間距離較遠(yuǎn)，實際上垂向速率陡變的地方大多在盆嶺地形的陡變帶上，而陡變帶難以布設(shè)測量點進(jìn)行控制，靠純數(shù)學(xué)的方法內(nèi)插很難給出理想的擬合結(jié)果?；谏鲜隹紤]，針對圖4的速率矢量數(shù)據(jù)進(jìn)行進(jìn)一步的融合計算，再給出垂直形變等值線。具體方案如下：

1)盆嶺之間垂直形變零值線的確定。根據(jù)點位已有速率分布和構(gòu)造地形分布，首先確定盆嶺構(gòu)造間垂直形變速率的陡變帶。因地形差異就是長期垂向構(gòu)造運動差異造成的，而水準(zhǔn)覆蓋時間都在10 a以上，這樣確定的垂直形變陡變帶比計算機在不考慮構(gòu)造地形的情況下均勻內(nèi)插的結(jié)果更合理、更接近實際情況。零值線就從陡變帶的中間穿過，將該零值線作為約束再繪制等值線圖，就會得到更接近實際情況的形變等值線圖。

2)在零值線的約束下，分別用最小二乘配置和多核面函數(shù)方法內(nèi)插擬合圖4的速度矢量數(shù)據(jù)。

3)最小二乘配置方法有利于提取低頻形變信息，特殊地形造成點位空間密度差異較大時其削峰填谷作用較強，容易失真；多面函數(shù)利用核點作控制，利于擬合陡變帶的高頻形變信息。為使二者互補，基于多核函數(shù)給出經(jīng)典平差分量(用于替代最小二乘配置中經(jīng)典平差部分)建立最小二乘配置的誤差方程，再利用赫爾默特迭代方法確定水準(zhǔn)和GNSS兩類資料的最佳權(quán)重，得到最佳的帶權(quán)平均結(jié)果作為最后的融合計算結(jié)果。

4)利用融合計算結(jié)果繪制垂直形變等值線圖(圖6，單位mm/a)。由圖可見，最終融合結(jié)果的垂直形變速率等值線圖與實際地形構(gòu)造分布的吻合度較高，較好地反映了該區(qū)域地殼垂向的繼承性構(gòu)造活動。

圖6 天山構(gòu)造區(qū)垂直形變速率等值線圖 (最終融合的結(jié)果，墨氏投影)Fig.6 Vertical deformation rate isoline map of the Tianshan tectonic region (final fusion results, Mercator projection)

5 討 論

因篇幅限制，本文就幾個問題進(jìn)行初步的討論。

1)GNSS與水準(zhǔn)垂向速率的差異。GNSS和精密水準(zhǔn)給出的高程值一個是法線上的高程，一個是垂線上的高程，二者的差異包括大地水準(zhǔn)面差距和垂線偏差對高程的影響。理論研究表明，大地水準(zhǔn)面差距雖未知，但對于同一個點在相當(dāng)長的時間內(nèi)可近似認(rèn)為是一個常數(shù)，其對點位垂向運動速率計算的影響可以忽略；垂線偏差通常是一個很小的角度，對垂向高程變化速率的影響也很微小[3]，可以忽略不計。

實際觀測和計算研究發(fā)現(xiàn)，在同一區(qū)域GNSS和水準(zhǔn)得到的垂向速率差異有時是比較大的。這一差異是多種因素的綜合影響，至于哪一個結(jié)果更可靠，還要具體問題具體分析。GNSS垂向坐標(biāo)測量精度遠(yuǎn)低于水平向坐標(biāo)的測量精度，其中氣象要素中的周期性干擾影響比較顯著，但多年觀測的平均垂向運動速率可較好地消除周期性干擾的影響；精密水準(zhǔn)的近距離兩點高差測量精度較高，對計算垂向形變速率在空間上的變化梯度有較好的控制作用，但高程和垂向運動速率觀測和計算結(jié)果中的誤差則會隨著水準(zhǔn)路線長度的增加而累積傳播。因此，對地殼垂直形變研究來說，二者各有千秋，綜合利用，進(jìn)行融合計算，可達(dá)到有益互補的作用。

2)同區(qū)域、同類資料獲得的空間形變場的差異。在同一區(qū)域，不同的學(xué)者使用同一類資料和相同的數(shù)據(jù)處理方法，也會得到不同的形變場圖像，讓人感到困惑。這個問題的根源是，資料類別雖然相同，但資料來源可能不同；甚至資料來源相同，但經(jīng)過預(yù)處理和篩選，最后采用的資料并不完全一樣。

3)相同資料不同數(shù)據(jù)處理方法獲得的形變場的差異。有時發(fā)現(xiàn)，所用資料完全一樣，但不同學(xué)者給出的形變場也存在一定的差異，進(jìn)一步比對發(fā)現(xiàn)，常與所用的數(shù)據(jù)處理方法、繪制形變場圖的軟件或采用的相關(guān)參數(shù)不同有關(guān)。比如測量平差的方法、程序、參數(shù)，數(shù)據(jù)內(nèi)插的方法、模型、軟件和參數(shù)等，每個環(huán)節(jié)有變化，最后得到的形變場圖都會存在一定的差異，有時差異還比較顯著。

4)地殼形變在地震預(yù)報中的意義。從前面的討論得知，不同的資料得到的形變場圖肯定不同；而資料相同，用不同方法、程序或參數(shù)得到的形變場圖像也會不同，因此，不同學(xué)者給出的形變場研究結(jié)果往往存在一定的差異。這會給用形變圖進(jìn)行地震預(yù)測研究帶來困惑，尤其是相關(guān)專業(yè)的研究人員在引用形變研究結(jié)果時，往往感到無所適從。

地殼結(jié)構(gòu)十分復(fù)雜，地殼應(yīng)力場發(fā)生變化時，因地殼介質(zhì)各向異性，又被很多尚未認(rèn)識到的斷層切割，地殼形變的響應(yīng)也必然十分復(fù)雜。從頻譜分析的角度講，地殼形變信息的頻譜成分十分豐富，不同測點的觀測資料包含的空間信息成分不同；不同觀測周期的資料包含的時間頻譜成分不同；不同的數(shù)據(jù)處理方法(各種平差、擬合、融合計算等)，其數(shù)學(xué)模型在某種程度上都能起到濾波的作用，對不同頻譜的形變信息有不同的選擇或抑制作用，有的側(cè)重于低頻長周期時空形變信息，有的側(cè)重于中頻信息，也有些側(cè)重于高頻信息，所以處理結(jié)果都不相同。哪一種數(shù)據(jù)處理方案最好，難以給出定論，因為不同構(gòu)造環(huán)境下孕育不同地震的形變異常信息在時間和空間上的頻譜成分不可能一樣。

因為測量點位的時空分布不一樣，所以研究者用不同資料、不同方法給出的地震形變異常圖像會存在差異，這是必然的結(jié)果。但基于同一空間分布的同一種資料，用同一套數(shù)據(jù)處理方法，在同一地區(qū)給出長時間的地形變研究的系列成果，往往會發(fā)現(xiàn)與地震對應(yīng)的異常變化，這對研究地震形變前兆具有十分重要的意義。根據(jù)多年的研究經(jīng)驗，類似本文給出的形變圖對中長期形變異常信息有比較好的體現(xiàn)，在幾年到幾十年的強震中長期預(yù)測研究中可起到較為重要的作用[7-14]。