橢球坐標(biāo)的多面函數(shù)法在區(qū)域地殼形變監(jiān)測中的應(yīng)用研究

胡 洋，王 璐

(西安市勘察測繪院，陜西西安710054)

一、引 言

隨著各個地區(qū)區(qū)域地殼形變監(jiān)測網(wǎng)的迅速發(fā)展，利用GPS觀測資料來研究區(qū)域地殼整體運動位移變化與內(nèi)部連續(xù)形變成為可能。但GPS的觀測結(jié)果在WGS-84橢球定義的參考橢球面上解算地殼應(yīng)變率，因此求解橢球坐標(biāo)下的應(yīng)變場對于區(qū)域地殼形變監(jiān)測的特征分析有著非常重要的意義［1］。

二、橢球坐標(biāo)下多面函數(shù)法的應(yīng)變場

1.橢球面上的應(yīng)變張量

根據(jù)彈性力學(xué)，應(yīng)變張量矩陣是一個對稱矩陣，在二維情況下由3個獨立分量構(gòu)成。橢球坐標(biāo)下的應(yīng)變張量矩陣是根據(jù)坐標(biāo)轉(zhuǎn)換的方法來進(jìn)行，基本思想是求出坐標(biāo)轉(zhuǎn)換矩陣。首先在形狀固定下的旋轉(zhuǎn)橢球上任意一點處建立一個正交曲線坐標(biāo)系，通過該坐標(biāo)系坐標(biāo)方向與以該橢球中心為坐標(biāo)原點的固定笛卡爾坐標(biāo)架的轉(zhuǎn)換關(guān)系，可推算出橢球面上活動坐標(biāo)系的位移梯度矩陣，從而根據(jù)該矩陣的線性組合與微分求解出旋轉(zhuǎn)橢球下的應(yīng)變張量［2］。因此可以將應(yīng)變張量矩陣的分量寫成式(1)的形式。

式中，h為大地高;λ為經(jīng)度;φ為緯度;Rn為卯酉圈曲率半徑;Rm為子午圈曲率半徑［2］。當(dāng)Rm=Rn=R，即橢球面變?yōu)榍蛎鏁r，(1)式就變?yōu)?/p>

式(2)就是應(yīng)變張量在球面坐標(biāo)系下解的形式。通過上述公式推導(dǎo)可知，應(yīng)變張量的球面解形式實際上是橢球面解形式的一種特定狀況［3］。由于地球的形狀可以近似看成一個橢球面，因而可以忽略式(1)中的大地高及高程變化信號，令uh=0、h=0整理得

2.多面函數(shù)模型

多面函數(shù)的基本模型如式(4)所示

式中，n為簡單數(shù)學(xué)面的張數(shù)或多層疊加面的層數(shù)，它的值與分塊擴充范圍內(nèi)參與點的個數(shù)相等;αj為待定參數(shù)，代表了第j個核函數(shù)對多層疊加面的貢獻(xiàn);R(x，y;xj，yj)為核函數(shù)。為了計算方便，多層疊加面中的n個核函數(shù)一般選用同一類型的簡單函數(shù)，通常是圍繞豎向軸旋轉(zhuǎn)的曲面，經(jīng)驗上在對地殼運動速率場的擬合中，選取下述函數(shù)作為核函數(shù)

式中，(x，y)為推估點的坐標(biāo);(xj，yj)為用于建立多面函數(shù)模型所選取的測點資料;δ為光滑因子;β為冪數(shù)，通過取不同的值來構(gòu)建不同的核函數(shù)，β一般取值為1/2(正雙曲面型)、3/2(三次曲面型)、-1/2(倒雙曲面型)。在建立多面函數(shù)模型時，對于研究區(qū)域內(nèi)已有的n個實測點數(shù)據(jù)。應(yīng)該在其中選取m個點數(shù)據(jù)(m≤n)來建立多面函數(shù)模型

式中，S為觀測值矩陣，因此有

3.橢球坐標(biāo)下多面函數(shù)應(yīng)變參數(shù)解算方法

對研究區(qū)域的GPS速度場的北向分量Vn與東向分量Ve，分別建立多面函數(shù)模型［3］

三、福建沿海地區(qū)應(yīng)變率場

為了進(jìn)一步研究福建地區(qū)地殼的水平形變特征分布，本文根據(jù)上述公式的推導(dǎo)進(jìn)行了Matlab編程實現(xiàn)，利用一組已有的GPS擬合速度場數(shù)據(jù)，采用0.2°×0.2°的經(jīng)緯格網(wǎng)密度，在扣除研究區(qū)域的整體剛性旋轉(zhuǎn)趨勢后得到其相應(yīng)的剩余速度分布場，以此來研究區(qū)域地殼內(nèi)部各處的相對形變特征。經(jīng)過交叉驗證法，建立出當(dāng) δ=1.14、β=3/2時，多面函數(shù)的建模結(jié)果最好。據(jù)此成功解算出基于多面函數(shù)法福建沿海應(yīng)變場，進(jìn)而由此來分析研究該區(qū)域現(xiàn)今地殼形變特征［4-7］。

如圖1所示，福建東南沿海地區(qū)最大、最小應(yīng)變率分布呈現(xiàn)出東南強于北西、東北強于南西的特征，其平均應(yīng)變值約為4.50×10-10/a;應(yīng)變率由沿海地區(qū)向內(nèi)陸延伸時由大變小，并且在斷裂帶交接區(qū)域值達(dá)到最大。結(jié)合該地區(qū)內(nèi)部相對形變場分布進(jìn)行分析，發(fā)現(xiàn)其特征十分相近。區(qū)內(nèi)主壓應(yīng)變優(yōu)勢方向較為一致，呈現(xiàn)出旋轉(zhuǎn)方向上的NWW-SEE。主張應(yīng)變軸走向近似水平，也呈現(xiàn)出順時針旋轉(zhuǎn)方向NE-NEE。從應(yīng)變值大小來看，主張應(yīng)變大于主壓應(yīng)變。在東南沿海地區(qū)閩江斷裂、莆田—南日島斷裂以及九龍江中下游斷裂與長樂—南澳斷裂帶的交匯地區(qū)，主張應(yīng)變均顯著增強，且成為該區(qū)域地殼的主應(yīng)變場。結(jié)合在NW與NE向斷裂的交匯處發(fā)育有福州盆地、漳州盆地及泉州盆地，且各盆地主軸線走向均為NW向，說明這一區(qū)的地形發(fā)育主要是由NW向的斷裂構(gòu)造帶進(jìn)行布控，且活動性應(yīng)該以拉張為主［8］。

圖中所示研究區(qū)域內(nèi)地殼主應(yīng)變率量級總體變化不大，只存在兩個較為明顯的突變區(qū):①福建西北部地殼水平形變量明顯減弱，其主應(yīng)變量級約為1.00×10-10/a 至 6.40×10-11/a。該區(qū)域的微弱應(yīng)變場范圍大，且過渡緩慢。結(jié)合前文的論述，可以推斷該區(qū)域現(xiàn)今地殼活動最為穩(wěn)定。②東南部長樂—南澳斷裂帶與北西向斷裂的交匯地區(qū)，地殼應(yīng)變率再次出現(xiàn)大幅度減弱趨勢，應(yīng)變率約為1.94×10-10/a至6.76×10-11/a，具有空間范圍小、衰減幅度大的特征;且這一衰減特征普遍適用于福建東南沿海長樂—南澳斷裂帶與主要的北西向斷裂交匯部位，尤其是在九龍江中下游斷裂的交匯部位，主壓應(yīng)變、張應(yīng)變在空間內(nèi)變化劇烈且急促，最大主壓應(yīng)變就分布在這一區(qū)域，其值約為4.43×10-8/a。應(yīng)變率場分布顯示出福建東南沿海區(qū)域地殼應(yīng)力場并不均勻，且其東南變化幅度大于北東。

圖2所示為福建地區(qū)地殼面膨脹值分布，其正值部分表示區(qū)域地殼處于膨脹拉張的狀態(tài)，相對應(yīng)的負(fù)值分布則表示該區(qū)域處于地殼的擠壓縮短之中。研究區(qū)域地殼引張最為劇烈的區(qū)域位于漳州西北側(cè)的九龍江斷裂與長樂—南澳斷裂帶交匯區(qū)域，約為2.90×10-8/a;而在長樂—南澳斷裂帶的最南段，地殼壓縮率達(dá)到最大值3.89×10-8/a?？傮w上看，福建地區(qū)壓縮區(qū)域分布呈帶狀顯示，近似北東平行于海岸線呈條帶狀連續(xù)分布;而在長樂—南澳斷裂與北西向斷裂的交匯區(qū)域依次出現(xiàn)相對離散化的局部引張，且地殼擴張的影響范圍幾乎錯斷了地殼縮短梯度帶，進(jìn)而形成了沿海岸線分布的地殼擴張—壓縮交替帶，在這一類地區(qū)地殼受力不均且變化較大。

圖1 最大、最小主應(yīng)變率

圖2 面膨脹

圖3為福建及其鄰近地區(qū)地殼南北向地殼應(yīng)變率，圖中以九龍江中下游斷裂帶為分界，地殼南北向應(yīng)變率出現(xiàn)顯著的差異性。斷裂以南，應(yīng)變梯度近似東西向呈條帶狀展布，地殼擠壓與拉張交替出現(xiàn)，且地殼形變趨勢明顯強于斷裂以北。在漳州以北地區(qū)，南北向應(yīng)變率梯度帶近似垂直于海岸線分布。由此可以反映出研究區(qū)域地殼形變南西強于北東，且該區(qū)域地殼應(yīng)力場并不均勻。東西向應(yīng)變率場如圖4所示，呈現(xiàn)出的特征與面膨脹頗為相似，說明造成區(qū)域地殼壓縮或膨脹差異的主要是來自近似南北向的應(yīng)力分量。結(jié)合GPS實測速度場與內(nèi)部相對形變場出現(xiàn)的差異可以推斷，福建及其鄰近地區(qū)地殼形變雖然以來自印度板塊北東向推擠產(chǎn)生的東向傳遞作用為主，但是其相應(yīng)應(yīng)力分量較為均勻。而來自菲律賓板塊的北西向推擠產(chǎn)生的應(yīng)力場則空間變化比較明顯，可能是該區(qū)域地震孕震的主要構(gòu)造因素。

圖3 南北向應(yīng)變率

圖4 東西向應(yīng)變率

圖5所示的研究區(qū)域最大剪應(yīng)變普遍較小，且近似呈NE-SW的條帶狀分布。這與長樂—南澳斷裂帶現(xiàn)今構(gòu)造活動呈現(xiàn)出的左旋壓扭性特征發(fā)育結(jié)果和該地區(qū)最大剪應(yīng)變較為一致。剪應(yīng)變最大值出現(xiàn)在長樂—南澳斷裂帶與九龍江斷裂帶的交匯地區(qū)，其值達(dá)到 8.99×10-8/a。

圖6為福建地區(qū)區(qū)域第一剪應(yīng)變率場分布，圖中正值區(qū)域表示該地區(qū)受到東西向伸長、南北向壓縮的形變;為負(fù)值時則相反。從圖中可以看出，福建地區(qū)從整個呈現(xiàn)由NW-SE方向的條段狀展布。西北部地殼比較穩(wěn)定，中部東西向伸長、南北向壓縮，沿海福州、漳州城市群東西向壓縮、南北向伸長，與該地區(qū)相對形變場一致。圖7為福建地區(qū)第二剪應(yīng)變分布圖，為正值時，表示該地區(qū)收到NE-SW向伸長、NW-SE向壓縮的形變;負(fù)值則相反。從圖中可推斷，造成福建地區(qū)的地殼壓縮或膨脹差異的還受到NE-SW向的應(yīng)力分量作用，在福州盆地地殼NESW向壓縮、NW-SE向伸長，南部的漳州盆地則相反，因而可以判斷福建地區(qū)地殼運動的受力比較復(fù)雜，差異性較高。

圖6 第一剪應(yīng)變

圖7 第二剪應(yīng)變

四、結(jié) 論

本文基于橢球坐標(biāo)下的多面函數(shù)法建立了福建沿海地區(qū)的多面函數(shù)模型，在0.2°×0.2°格網(wǎng)節(jié)點GPS擬合速率場上，通過交叉驗證法，認(rèn)為當(dāng)光滑因子δ=1.14、核函數(shù)為β=3/2時，建模效果最優(yōu)。根據(jù)此多面函數(shù)模型進(jìn)一步導(dǎo)出了基于旋轉(zhuǎn)橢球面上的應(yīng)變張量公式，并解算得到了福建沿海地區(qū)的各種應(yīng)變場。

以往對于福建及其鄰近地區(qū)沿海地殼形變特征的研究成果認(rèn)為:該區(qū)域地殼形變存在東強西弱，南強北弱的總體特征。根據(jù)本文的研究結(jié)果，研究區(qū)域地殼的現(xiàn)今水平形變特征可能更加復(fù)雜。在總體形變的背景場相一致的基礎(chǔ)上，局部存在較大的差異性變化。具體地，以助漳州盆地及九龍江中下游斷裂為分界，可以將研究區(qū)域地殼水平形變分為兩個小區(qū)。在九龍江斷裂以北的沿海地區(qū)，地殼水平形變特征較為一致，基本符合南強北弱、東強西弱的特點。根據(jù)以上各個應(yīng)變場的特征分析，為進(jìn)一步研究區(qū)域地殼形變監(jiān)測的深部驅(qū)動機制提供參考。

［1］黃聲享，尹暉，蔣征.變形監(jiān)測數(shù)據(jù)處理［M］.武漢:武漢大學(xué)出版社，2003.

［2］劉序儼，黃聲明，林巖釗.地形旋轉(zhuǎn)張量探討［J］.大地測量與地球動力學(xué)，2010，30(5):57-63.

［3］劉序儼，黃聲明，梁全強.旋轉(zhuǎn)橢球面上的應(yīng)變與轉(zhuǎn)動張量表達(dá)［J］.地震學(xué)報，2007，29(3):240-249.

［4］吳紹祖，丁學(xué)仁，李祖寧，等.福建地區(qū)現(xiàn)今地殼運動變化特征研究［J］.華南地震，2005，25(3):87-94.

［5］石耀林，朱守彪.用GPS位移資料計算應(yīng)變的方法的討論［J］.大地測量學(xué)與地球動力學(xué)，2006，26(1):1-8.

［6］劉經(jīng)南，施闖，姚宜斌，等.多面函數(shù)擬合法及其在建立中國地殼平面運動速度場模型中的應(yīng)用研究［J］.武漢大學(xué)學(xué)報:信息科學(xué)版，2001，26(6):501-503.

［7］張貴鋼，楊志強，王慶良，等.基于最小二乘配置模型探討青藏高原巴彥喀拉塊體的運動特征［J］.地球科學(xué)與環(huán)境學(xué)報，2011，33(2):213-216.

［8］CHEN Z，BCB，LIU Y，etal.Global Positioning for India/Eurasia intercontinental deformation［J］.Geophys Res，2001(105):16215-16227.