楊瑩輝，陳強＊，劉國祥，程海琴，2，劉麗瑤，胡植慶

1西南交通大學遙感信息工程系，成都 610031

2華東交通大學土木建筑學院，南昌 330013

3臺灣大學地質(zhì)科學系，臺灣臺北

1 引言

2008年5月12日，以汶川縣映秀鎮(zhèn)為震中發(fā)生了MS8.0級特大地震，地震發(fā)生機制為龍門山中央斷裂的右旋走滑逆沖，其動力誘因為印度板塊持續(xù)向北運動擠壓歐亞大陸板塊，使得青藏高原東緣地塊向東運動，而與之相鄰的四川盆地地塊較為堅固（徐錫偉等，2005；單新建等，2009），導致龍門山斷層應力持續(xù)積累并趨于臨界狀態(tài)，最終積蓄能量從映秀—北川主斷裂帶釋放導致汶川地震的發(fā)生.

汶川地震引發(fā)龍門山地區(qū)強烈的地表形變，同震形變場的準確獲取對于分析龍門山地區(qū)地表破裂的空間分布特征、震源參數(shù)計算以及地球物理反演具有重要的科學意義.汶川地震后，中國地殼運動觀測網(wǎng)絡項目組及時對龍門山地區(qū)的GPS觀測站進行了復測，并公布了觀測結果（張培震，2008）.GPS形變觀測顯示，以映秀—北川斷裂帶為中心，兩側斷層上下盤存在顯著的相向運動和水平縮短現(xiàn)象，同時斷層下盤以下降為主，上盤在地表破裂線周邊呈現(xiàn)上升趨勢，但在遠離破裂線后很快又轉為下降運動，GPS形變的這些特征揭示了斷層運動以逆沖為主的特性（張培震，2008）.但由于GPS稀疏分布的特性，難以獲得發(fā)震區(qū)完整的高分辨率地表形變場，不利于斷層面破裂的精細反演計算.汶川地震前后，ALOS衛(wèi)星搭載的PALSAR傳感器獲取了該區(qū)域大范圍SAR影像數(shù)據(jù)，InSAR技術可為龍門山地區(qū)提供大范圍、高分辨率的地表形變觀測資料（單新建等，2009）.

汶川地震InSAR同震形變場的觀測研究表明，由于雷達衛(wèi)星數(shù)據(jù)受軌道誤差、大氣信號延遲等負面因素的影響，從InSAR干涉圖中提取精細化的地表形變場干擾較大，尤其是汶川地震波及范圍廣，鄰軌PALSAR干涉圖之間存在嚴重的形變場不連續(xù)問題（Zhang etal.，2008；Hashimoto etal.，2009；Shen etal.，2009；Furuya etal.，2010；Liu etal.，2010；Tong etal.，2010），積極探索鄰軌InSAR同震位移的修正方法，有助于準確獲取龍門山地區(qū)精細的同震形變場與斷層滑移參數(shù)的精化反演.

為解決InSAR軌道誤差殘留相位與鄰軌形變條紋的不連續(xù)問題，本文提出聯(lián)合GPS形變觀測與鄰軌相位修正的數(shù)據(jù)處理策略，利用稀疏分布的GPS位移觀測去除PALSAR軌道誤差引入的殘留平地相位，采用鄰軌約束平滑消除干涉圖之間的不連續(xù)問題，從而獲得龍門山地區(qū)大范圍精細的同震形變場，在此基礎上，采用彈性半空間位錯模型反演地殼發(fā)震斷層的滑移模型，并由此計算汶川強震的地震矩數(shù)值.

2 ALOS／PALSAR影像干涉處理

使用ALOS衛(wèi)星PALSAR傳感器獲取的地震前后雷達影像開展干涉處理，以其中連續(xù)6個升軌條帶共46幅影像作為同震形變觀測數(shù)據(jù)，各軌道影像的基本信息如表1所示.

雷達影像各條帶的場景覆蓋范圍如圖1紅色矩形方框所示，圖中的綠線為野外調(diào)查獲得的汶川地震地表破裂線，主要包括中央主斷裂（映秀—北川）、小魚洞斷裂和前山斷裂（灌縣—彭縣），小魚洞斷裂帶位于中央斷裂帶和前山斷裂帶之間，前山斷裂帶長60余公里，在地下約10km處并入中央主斷裂帶（李勇等，2009）.圖1中的藍色與棕色箭頭表示龍門山地區(qū)的GPS同震形變觀測向量（張培震，2008），其中棕色箭頭代表單位長度形變?yōu)?0cm，藍色箭頭單位長度形變?yōu)?0cm.圖中的紅色小圓為MS＞4.0級的余震震中位置（Zhang etal.，2008；Hashimoto etal.，2009；Shen etal.，2009；Furuya etal.，2010；Liu etal.，2010；Tong etal.，2010）.

表1 汶川地震PALSAR影像的干涉對信息Table 1 PALSAR interferometric pairs used in the study of Wenchuan earthquake

對上述PALSAR影像進行干涉處理，使用ASTER DEM作為外部DEM以去除地形相位分量，獲得雷達差分干涉相位圖（如圖2）.從圖2可以看出，相鄰軌道間的干涉條紋存在顯著的不連續(xù)性以及軌道殘留相位的干擾，在475—476軌道鄰接邊界（左上角的黑色方框）存在明顯的干涉條紋不連續(xù)，而圖幅中間黑色方框所標識的473條帶則存在顯著的軌道殘留相位，導致遠場存在明顯的干涉條紋.此外，在右下角的471—472條帶連接處，則存在嚴重的大氣相位，其形變條紋與473—474存在顯著的不一致性，造成形變場的不連續(xù).為克服上述干涉條紋形變場的不連續(xù)問題，本文提出采用GPS形變位移對PALSAR軌道誤差殘留平地相位進行建模修正，利用鄰軌平滑約束條件校正鄰軌干涉條紋.

3 軌道殘留相位建模與校正

在InSAR數(shù)據(jù)處理過程中，平地相位（也稱為參考相位）的計算是基于SAR軌道狀態(tài)矢量（包括位置矢量與速度矢量），而ALOS衛(wèi)星軌道數(shù)據(jù)的不確定性導致干涉圖中存在殘留平地相位，在圖2中473條帶的干涉相位表現(xiàn)尤為顯著，在473條帶的黑色矩形框內(nèi)，地震遠場呈現(xiàn)出顯著的周期性條紋，從形變場空間分布特征表明，該信號為殘留的平地相位.

圖1 汶川地震PALSAR影像場景覆蓋范圍Fig.1 Ground coverage of PALSAR images for the Wenchuan earthquake

圖2 InSAR同震形變干涉相位Fig.2 InSAR derived interferogram of coseismic deformation field

圖3 軌道殘留相位校正前后的干涉相位圖（a）校正前的干涉相位；（b）校正后的干涉相位.Fig.3 Local interferogram before（a）and after（b）correction from residual orbit phase modeling

為消除軌道數(shù)據(jù)誤差引入的殘留參考相位，本文提出以GPS精密形變測量值作為約束，對同一區(qū)域的InSAR與GPS形變觀測進行差分建模，從而消除或減弱由于衛(wèi)星軌道誤差引入的殘留相位值.

首先將GPS觀測得到的三維形變轉換至雷達視線方向上的LOS形變（陳強等，2012a），定量的數(shù)學轉換模型如下：

式中，gLos代表GPS三維形變投影到衛(wèi)星LOS視線方向的形變，dUp、dN、dE分別為GPS測量得到的豎直向、南北向和東西向形變量，θ代表衛(wèi)星側視角，α是衛(wèi)星飛行方向與北方向的夾角，θ和α可依據(jù)衛(wèi)星參數(shù)予以計算.

從圖2中未校正的473條帶可以看出，軌道殘留相位在空間尺度上呈現(xiàn)出顯著的線性分布（張磊等，2007），因此構建如下的軌道殘留誤差修正函數(shù)模型：

式中，gLos為GPS形變轉換到LOS方向的形變，可通過式（1）計算得到，dLos為InSAR觀測得到的形變，l、p分別表示GPS測站在雷達干涉圖中的方位向（Azimuth）和距離向（Range）坐標，a0、a1、a2為待求模型系數(shù)，采用3個或3個以上的GPS形變觀測量，即可求解出上述系數(shù)，從而建立軌道殘留相位修正模型.

從圖2中的三維GPS形變測站分布來看，PALSAR條帶的472、473、474、475四個條帶中，GPS點數(shù)量均大于3個，可用于軌道殘留相位的建模，而472條帶中的H010號GPS點位于余震分布密集的區(qū)域，其InSAR相干性較差，形變觀測數(shù)據(jù)質(zhì)量不高，不用于軌道誤差的校正.

采用上述模型和GPS觀測數(shù)據(jù)分別對473、474和475條帶中的殘留相位予以建模，獲得校正后的差分干涉圖，如圖3所示（對應于圖2中的矩形框區(qū)域），對比校正前后干涉相位，可以看出由于軌道不確定性所導致的殘留相位在校正后相位圖中已得到較好的消除，在473—474條帶的邊界已經(jīng)沒有明顯的多余條紋.

4 鄰軌平滑校正

對于PALSAR條帶內(nèi)無足夠數(shù)量GPS校正點的情況下，如475—476的條帶連接處，本文提出了基于鄰軌間干涉形變差異建模的相對校正方法.由于主震震級造成的地表形變遠大于余震和其他因素引起的形變（如震后余滑和孔隙回彈效應在短時間內(nèi)造成的形變量較小，與主震形變量級相比不顯著），因此可以認為在上述PALSAR影像的時間基線跨度內(nèi)，相鄰條帶重疊數(shù)據(jù)具有平滑一致的形變量.以GPS校正后的干涉圖為基準，對鄰軌干涉圖與基準干涉圖之間的形變差異進行建模，從兩幅鄰軌重疊影像中提取高相關性的目標點，解算形變差異模型系數(shù)，進而實現(xiàn)對未校正PALSAR條帶形變場與干涉相位的修正.該方法充分利用鄰軌重疊區(qū)域已校正的干涉相位數(shù)據(jù)，對無足夠可用數(shù)量GPS形變的干涉條帶進行校正，可進一步減弱未校正條帶內(nèi)的軌道相位殘差和大氣殘留相位.

圖4為鄰軌條帶的振幅影像，藍色虛線方框為相鄰條帶間的重疊部分，顧及條帶重疊區(qū)的同震形變應保持一致性的原則，以其中一幅經(jīng)GPS校正后的干涉數(shù)據(jù)作為基準，對另一幅未校正干涉條帶相位值進行建模修正.例如，475條帶已經(jīng)過GPS校正，且與474和473形變條紋保持較好的一致性，則以475為基準對鄰軌的476干涉條紋進行建模修正，具體算法如下.

首先在鄰軌PALSAR影像的重疊區(qū)進行同名點識別，基于相關系數(shù)最大值匹配方法開展重疊區(qū)同名點搜索，結果如圖4所示，圖中的藍色和白色點代表所識別的同名點.為保證后續(xù)建模的可靠性，提取相關系數(shù)大于0.4的同名點（圖中白色點）用于數(shù)學建模，計算高相關性同名點在兩個條帶中的形變差，依據(jù)其在未校正條帶中的方位向和距離向坐標，構建如下式的形變校正模型：式中，Δd表示兩條帶同名點形變量之差，l、p分別為同名點在被校正條帶中的方位向和距離向坐標，bi（i＝0，1，2，3，4，5）為待求模型參數(shù).當高質(zhì)量的同名點數(shù)量大于6，使用最小二乘方法解算出模型系數(shù)，然后使用該模型對鄰軌條紋相位進行數(shù)據(jù)校正，獲得校正后的差分干涉圖，圖5為條帶475與476經(jīng)鄰軌校正后的干涉相位圖，從圖5中可以看出，條帶之間的不連續(xù)性問題已得到較好解決，其鄰軌干涉條紋連續(xù)性較好.

對上述校正后的差分干涉圖進行相鄰條帶相位數(shù)據(jù)拼接，獲得完整的地表形變干涉相位圖，如圖6所示，使用Snaphu軟件對拼接后的干涉相位進行全局相位解纏和絕對形變解算（Chen and Zebker，2001；Chen and Zebker，2002；陳強等，2012b），獲得地表沿衛(wèi)星雷達視線方向（LOS）的形變場，顧及圖2中的JYAN GPS基站所在的InSAR像元具有較高的相干性，且遠離斷層破裂線，選取該GPS觀測點作為相位解纏的形變基準點，最終獲得PALSAR影像區(qū)域完整的同震形變場，如圖7所示.

圖4 鄰軌PALSAR振幅影像重疊區(qū)與聯(lián)系點Fig.4 Spatial overlapping of adjacent PALSAR amplitude maps and detected tie points

圖5 鄰軌約束校正前后的干涉相位（a）校正前的干涉圖；（b）校正后的干涉圖.Fig.5 Local interferogram before（a）and after（b）correction using adjacent geometric constraints

圖6 —7所示的校正后InSAR同震形變場與已有研究成果相比（單新建等，2009；Shen etal.，2009；Furuya etal.，2010；Liu etal.，2010），整體趨勢基本一致，例如，觀測到最大的地表形變場位于北川區(qū)域，LOS向形變量級在1.2m左右，汶川—映秀段斷層上盤存在量級0.6～0.8m的顯著形變區(qū)，成都平原存在厘米量級的地表形變（Hashimoto etal.，2009；Liu etal.，2010；Shen etal.，2009；Furuya etal.，2010；Tong etal.，2010）.

然而，校正后的InSAR形變場在局部細節(jié)區(qū)域呈現(xiàn)出顯著的不同，首先，已有研究中條帶473存在的嚴重軌道誤差在本文中通過GPS校正方法得以很好去除，并在此基礎上，使用473條帶的觀測數(shù)據(jù)對條帶472和471進行了鄰軌平滑校正，顯著減弱了條帶471—472中存在的大氣、軌道等系統(tǒng)誤差.此外，依據(jù)鄰軌校正方法，條帶476也基于475進行了誤差校正.校正后的干涉相位較好地解決了已有研究成果中顯著存在的鄰軌拼接不連續(xù)現(xiàn)象，使得最終提取的形變場具有較好的鄰軌平滑性，有利于后續(xù)斷層滑動分布模型的精確反演.

圖6 校正后的InSAR同震形變干涉條紋Fig.6 Corrected InSAR interferogram related to coseismic deformation

圖7 校正后的汶川地震InSAR同震形變場Fig.7 Corrected InSAR surface deformation field of Wenchuan earthquake

5 同震形變場精度評定

5.1 絕對精度分析

為定量分析上述校正后InSAR同震形變場的精度與可靠性，以高精度GPS觀測值作為參考，對InSAR形變與同點位GPS形變進行差異比較，統(tǒng)計分析校正后InSAR同震形變場的實際絕對精度.

將PALSAR干涉圖中的GPS點位分為兩類，干涉相干性在0.3以上的點位劃分為高相干點，低于該相干閾值則為低相干點，將GPS觀測得到形變位移投影轉換到雷達LOS方向上的位移，定量數(shù)據(jù)比較如表2所示.

從表2中可以看出，經(jīng)軌道殘留相位校正和鄰軌平滑校正后的同震形變場，高相干點的形變量殘差由8.95cm減小至3.59cm，精度提高了60%，低相干點平均殘差由25.21cm減小至15.00cm，精度提高達40%，表明本文提出的軌道殘留相位建模和鄰軌平滑方法具有較高的可靠性，同震形變場的精度得以顯著提高.

5.2 鄰軌形變的相對精度分析

對于地震斷層錯動引發(fā)的地表形變，在斷層上下盤錯動的邊界處通常存在顯著的形變突變，如逆斷層上盤呈抬升運動，而下盤則表現(xiàn)為下降，對于其余方位的形變場則較為平滑.因此，在InSAR同一個條帶內(nèi)的同震形變場，除了斷層上下盤分界處的地表破裂帶區(qū)域形變不連續(xù)之外，其余地方則呈現(xiàn)出較好的平滑度和連續(xù)性.然而，在相鄰兩個不同PALSAR條帶之間的邊界區(qū)域，則可能出現(xiàn)前述的由于軌道誤差、大氣相位等因素引起的干涉相位或形變場不連續(xù)現(xiàn)象，因此采用條帶邊界形變的平滑程度指標，可對同震形變場的可靠性與平滑度進行相對精度評價.

本文引入拉普拉斯平滑參數(shù)作為同震形變場連續(xù)性的定量評價指標，以中心像元S（i，j）為計算單元，通過鄰近像元的差分計算，獲得干涉形變場的拉普拉斯二階導數(shù)，如下式所示：（4）式中，i，j為待計算像素位置，S（i，j）為像元位置處的形變量，Δl1和 Δl2為SAR分辨單元尺寸，ρi，j為平滑因子，｜ρi，j｜的數(shù)值越小，表示形變場的平滑度越高（Jonsson etal.，2002）.

表2 修正前后的InSAR形變與GPS形變的差異比較Table 2 Displacement comparison between InSAR and GPS derived deformation field

對校正前后InSAR形變場各相鄰條帶邊界處的相位，使用上述公式（4）計算其平滑度指標，得到如圖8所示的鄰軌邊界平滑指數(shù)分布.

從圖8可以看出，經(jīng)軌道誤差和鄰軌約束校正后，條帶邊界干涉形變場的非零｜ρi，j｜像元數(shù)量顯著減少，邊界平滑度得到明顯改善，對所有鄰軌邊界像元平滑因子ρi，j進行數(shù)據(jù)統(tǒng)計，直方圖如圖9所示，統(tǒng)計數(shù)據(jù)表明，校正后形變場平滑因子多趨于零，標準偏差減小了33%.對平滑因子進行分區(qū)間統(tǒng)計，在平滑因子絕對值較小的（-40，40）區(qū)間內(nèi)，校正后形變場平滑因子比校正前提高了7.7%，對于平滑因子絕對值較大的（-200，-40）、（40，200）區(qū)間，校正后形變場的平滑因子減少了6.8%.綜合以上分析，經(jīng)過軌道誤差校正和鄰軌約束校正后的形變場，鄰軌邊界具有更好的平滑度，形變條紋的不連續(xù)性得以較好修正.

圖8 同震形變場校正前后的鄰軌邊界平滑度（a）校正前的邊界平滑度；（b）校正后的邊界平滑度.Fig.8 The smoothness of coseismic interferogram respectively before（a）and after（b）correction

圖9 形變場的鄰軌邊界平滑因子統(tǒng)計（a）校正前鄰軌邊界平滑因子；（b）校正后鄰軌邊界平滑因子.Fig.9 Statistics of adjacent track smoothing factor of derived deformation field

6 斷層滑移模型反演

在精確獲得同震地表形變場基礎上，采用Okada彈性半空間位錯模型反演發(fā)震斷層表面的滑移分布（Okada，1985）.通過對野外地質(zhì)調(diào)查數(shù)據(jù)分析和綜合已有研究成果（Zhang etal.，2008；Hashimoto etal.，2009；Liu etal.，2009；Shen etal.，2009；Furuya etal.，2010；Liu etal.，2010；Tong etal.，2010），本文構建了如圖10中所示的F1—F6六段斷層分段模型，首先假設斷層表面為均勻滑動，使用非線性遺傳算法對斷層參數(shù)進行優(yōu)化反演，為保證反演結果的可靠性，我們賦予反演模型不同的初始值，進行多次反演搜優(yōu)計算，并對反演參數(shù)集進行統(tǒng)計分析，最終獲得了各子斷層的幾何參數(shù)及其對應精度，如表3所示.

表3 汶川地震斷層幾何Table 3 Fault geometry of Wenchuan earthquake

對于中央斷裂帶（映秀—北川）長大斷層來說，斷層面的實際滑動往往呈現(xiàn)為非均勻分布，為更精細獲得局部斷層滑動的空間特征，對表3中各子斷層進行4km×4km的格網(wǎng)離散化處理，共獲得1051個子斷層，以其作為基本單元精細反演斷層滑動的空間分布，采用拉普拉斯平滑約束斷層滑動模型，如式（5）所示：

式中，Gi為根據(jù)離散子斷層單位走滑量和逆沖量計算得到的InSAR數(shù)據(jù)反演格林函數(shù)，Hi為GPS數(shù)據(jù)反演格林函數(shù)，α、β分別對應其權重因子，Los、def為InSAR和GPS觀測得到的地表形變量，Li為拉普拉斯約束格林函數(shù)，κ代表平滑因子分別代表斷層的走滑量和逆沖量.

由于GPS形變觀測在垂直方向的精度較弱，且垂直向觀測資料較少，因此在反演過程中僅使用了GPS平面觀測數(shù)據(jù).InSAR形變觀測數(shù)據(jù)的條帶471和472在校正后仍可能存在電離層及大氣負面影響，為避免其對反演結果造成不利影響，471—472條帶中形變場連續(xù)性較差的數(shù)據(jù)不參與反演計算.此外，為盡可能減弱不可靠的形變觀測結果對反演的影響，使用InSAR相干性大于0.3的形變觀測數(shù)據(jù)，并對其進行四叉樹降采樣工作，該采樣方法在高形變梯度區(qū)域采樣量大，在低形變梯度區(qū)域采樣數(shù)據(jù)量小，能最大限度地保留形變場的空間特征.降采樣過程一方面可避免地震遠場InSAR非構造形變噪聲對反演結果的影響，另一方面也可有效避免因觀測數(shù)據(jù)過多引起的計算負荷過重和計算效率低下的問題，最終的降采樣形變場如圖10所示.

圖10 四叉樹降采樣InSAR形變場Fig.10 Distribution of InSAR sampled points based on quad-tree method

對于InSAR和GPS形變觀測數(shù)據(jù)用于斷層滑移反演的權重分配，本文將依據(jù)其觀測精度進行確定.根據(jù)前述的形變場精度評定可知，高質(zhì)量InSAR觀測數(shù)據(jù)的殘差為0.036m，在實際反演中使用相干性較高的觀測點，此處適當放寬其精度取為0.04m，對于GPS形變觀測數(shù)據(jù)的精度，根據(jù)相關研究可知其平面精度約為2～5mm，取3mm作為精度指標（張培震，2008）.因此，InSAR和GPS點的誤差之比為13.3∶1.0，反比例計算可得其相應權重分別為α＝0.075、β＝1.0.

在確定了數(shù)據(jù)權重之后，賦予反演模型連續(xù)變化的平滑因子κ，并依據(jù)反演結果計算獲得相應斷層面粗糙度與模型殘差的對應關系，如圖11所示，從圖中可看出，隨著粗糙度的增大，模型殘差急劇減小，顧及模型殘差和斷層粗糙度的折中關系，我們選取了如圖11中箭頭標明處的粗糙度為7.5cm／km處的平滑因子為最優(yōu)的κ取值.

圖11 模型殘差與斷層面粗糙度關系曲線Fig.11 The trade-off curve with misfit function of InSAR and GPS plotted as function of average solution roughness

圖12 斷層反演滑動及其精度分布（a）斷層滑動三維分布；（b）斷層滑動誤差分布；（c）斷層滑動平面分布.Fig.12 Fault slip distribution and its error（a）Three-dimensional distribution of coseismic slip；（b）Errors of coseismic slip；（c）Coseismic slip distribution on fault surface.

確定了各項模型參數(shù)后，在非負最小二乘原則下進行斷層滑動量的反演計算，獲得汶川地震斷層滑動的空間分布，如圖12所示，其中圖12a為反演斷層滑動的三維分布.為了驗證滑動反演結果的可靠性，我們基于蒙特卡洛算法生成100組帶隨機擾動的觀測數(shù)據(jù)集，并基于這些數(shù)據(jù)反演相應的斷層滑動分布，對所有計算得到的斷層滑動分布集進行統(tǒng)計分析，從而估算出反演模型的精度分布，如圖12b所示.統(tǒng)計表明反演模型的平均精度為0.014m，表明圖12a的反演結果是可靠的.同時需要指出的是，最大誤差分布位于青川地區(qū)16km深度處，分析此現(xiàn)象的可能原因，一方面是該區(qū)域的InSAR觀測數(shù)據(jù)由于存在大氣影響而未用于參與滑動反演，導致反演模型對此部分斷層滑動的約束較弱；另一方面是汶川主震后存在密集余震，InSAR觀測的滯后性使得相鄰條帶形變場中包含了不同余震的形變信息，對反演精度造成了一定影響.

為清晰顯示斷層面的滑動分布狀態(tài)，我們繪制了如圖12c所示的斷層滑動平面分布圖，從中可看出，斷層滑動整體表現(xiàn)為逆沖為主兼具右旋走滑特性，平面滑動量主要分布于映秀、北川和青川區(qū)域，斷層最大滑動量為9.0m，位于北川地表附近，導致北川地區(qū)地表破裂最為嚴重.在地殼深度方向，地表以下0～16km范圍內(nèi)孕育了較大的滑動成分，而16km以下滑動量分布并不顯著.

由于斷層滑動反演所使用的地表InSAR形變場約束以及斷層幾何模型等的不同，本文的反演結果與已有研究成果存在一定差異.首先，斷層北川段的滑動量在沿斷層走向分布更廣，但在深度方向上相對較淺，而映秀—汶川段的平均滑動量級較已有研究成果偏?。℉ashimoto etal.，2009；Shen etal.，2009；Furuya etal.，2010；Liu etal.，2010；Tong etal.，2010）；此外，青川段滑動量分布更為深入地下，這可能與青川段震后密集的余震分布相關（朱艾斕等，2008）.

根據(jù)上述反演出的斷層滑移模型與斷層幾何參數(shù)，正演計算汶川地震的InSAR地表形變場，如圖13中的左下插圖（a），與圖6中校正后InSAR同震形變場相比，兩者具有較高的一致性，尤其是在映秀和北川等形變嚴重的區(qū)域，正演數(shù)據(jù)與InSAR實際觀測值吻合的較好.但需要說明的是，在主斷層北端的青川區(qū)域，由于InSAR觀測數(shù)據(jù)存在嚴重大氣干擾負面影響，形變觀測值可靠性相對較弱，未參與到滑動分布模型的反演，使該部分數(shù)據(jù)與正演形變場存在顯著差異.

進一步分析同震形變場的模型反演殘差與精度，圖13為GPS的反演殘差向量（各矢量箭頭）與InSAR的反演殘差數(shù)值分布.圖（b）與（c）為GPS觀測值與模型計算值的對比圖，GPS在北東方向的反演殘差分別為4.3cm和6.2cm.此外，統(tǒng)計表明約有75%的GPS平面殘差小于2cm，且僅在近斷層附近存在殘差較大的現(xiàn)象.GPS的整體形變殘差為5.5cm，與同類已有研究基本相當.然而，本文InSAR反演形變殘差有顯著改善，從圖13中的InSAR殘差分布中可看出，大部分區(qū)域殘差小于10cm，形變殘差統(tǒng)計值達到9.2cm，與已有同類研究成果（平均殘差約15cm）對比（Zhang etal.，2008；Hashimoto etal.，2009；Liu etal.，2009；Shen etal.，2009；Furuya etal.，2010；Tong etal.，2010），InSAR形變反演精度提高約30%，驗證了InSAR同震形變場的高精度特性與可靠性，同時也表明前述InSAR形變場校正方法的有效性與準確性.

需要指出的是，斷層滑動模型的反演殘差一方面來源于地表形變場約束的觀測誤差（包括InSAR形變場中殘留的大氣、軌道誤差等），另一方面，主要受限于不確定的復雜斷層幾何、位錯理論模型的系統(tǒng)誤差等.本文斷層滑動反演的統(tǒng)計殘差為9.2cm，仍大于InSAR高質(zhì)量點3.6cm的觀測校正殘差，這主要是由于斷層幾何參數(shù)的精度、以及Okada位錯模型未考慮地殼不均勻特性所帶來的系統(tǒng)誤差.其中，斷層幾何對反演結果影響最為顯著，已有研究均假定斷層幾何在地殼不同深度具有相同的斷層傾角，這與實際斷層地下延伸模式可能存在差異.因而，即使在獲得高精度地表形變場的情況下，對于龍門山斷裂帶復雜的斷層幾何來說，仍可能出現(xiàn)相對較大的反演殘差.但相對于已有研究平均15cm的反演殘差，本文反演精度（9.2cm）有較大提高，這也正是InSAR形變場校正對于改善滑動模型反演精度的重要貢獻，所獲得的斷層滑動精細分布有助于進一步推進對汶川地震斷層破裂機制的認識和深入理解.

依據(jù)上述反演所獲得的斷層幾何參數(shù)和滑動分布模型，計算得到本次地震的地震矩為8.0469×1020N·m，對應的矩震級Mw7.9，與已有反演結果基本相當（Tong etal.，2010），但本文計算得到的地震矩稍顯偏大，一方面是因為本文使用了經(jīng)鄰軌校正后的精細化地表同震形變場，另一方面因為本文所用的斷層模型幾何傾角與已有方法存在一定差異.

圖13 InSAR和GPS形變的聯(lián)合反演殘差（a）模擬InSAR干涉圖；（b）GPS北方向形變與模型正演對比；（c）GPS東方向形變與模型正演對比.Fig.13 Residuals of InSAR and GPS joint inversion（a）Simulated InSAR interferogram；（b）Model prediction versus GPS observation in north direction；（c）Model prediction versus GPS observation in east direction.

7 結論

本文采用軌道誤差校正方法與鄰軌平滑技術對汶川地震形變場進行修正處理，獲得汶川地震精細的同震形變場.GPS形變觀測與InSAR同震形變的殘差分析表明，校正后形變場高質(zhì)量數(shù)據(jù)殘差為3.6cm，精度提高約60%，低質(zhì)量數(shù)據(jù)精度提高約40%；鄰軌邊界形變場平滑度分析表明，校正后形變場邊界平滑因子的標準差減小33%，驗證了本文提出的同震形變場校正方法的有效性和可靠性.

基于Okada彈性半空間位錯模型的斷層滑動反演結果表明，斷層滑動量在沿斷層走向上主要分布于映秀、北川和青川區(qū)域，在深度方向上主要集中在地殼淺部的0～16km范圍內(nèi)，最大滑動量為9.0m，位于北川縣城地表附近，InSAR反演后的模型殘差達9.2cm，反演精度約有30%的提高，由此計算出的地震矩數(shù)值為8.0469×1020N·m.

Chen C W，Zebker H A.2001.Two-dimensional phase unwrapping with use of statistical models for cost functions in nonlinear optimization.Journal of the Optical Society of America A，18（2）：338-351，doi：10.1364／JOSAA.18.000338.

Chen C W，Zebker H A.2002.Phase unwrapping for large SAR interferograms：statistical segmentation and generalized network models.IEEE Transactions on Geoscience and Remote Sensing，40（8）：1709-1719，doi：10.1109／TGRS.2002.802453.

Chen Q，Liu G X，Hu Z Q，etal.2012a.Mapping ground 3-D displacement with GPS and PS-InSAR networking in the Pingtung area，southwestern Taiwan，China.Chinese Journal of Geophysics（in Chinese），55（10）：3248-3258，doi：10.6038／j.issn.0001-5733.2012.10.007.

Chen Q，Yang Y H，Liu G X，etal.2012b.InSAR phase unwrapping using least squares method with integer ambiguity resolution and edge detection.Acta Geodaetica et Cartographica Sinica（in Chinese），41（3）：441-448.

Furuya M，Kobayashi T，Takada Y，etal.2010.Fault source modeling of the 2008Wenchuan Earthquake based on ALOS／PALSAR data.Bulletin of the Seismological Society of America，100（5B）：2750-2766，doi：10.1785／0120090242.

Hashimoto M，Enomoto M，F(xiàn)ukushima Y，etal.2009.Coseismic deformation from the 2008Wenchuan，China，Earthquake derived from ALOS／PALSAR images.Tectonophysics，491（1-4）：59-71，doi：10.1016／j.tecto.2009.08.034.

Jonsson S，Zebker H A，Segall P，etal.2002.Fault slip distribution of the 1999 Mw7.1Hector Mine，California，Earthquake，estimated from satellite radar and GPS measurements.Bulletin of the Seismological Society of America，92（4）：1377-1389，doi：10.1785／0120000922.

Li Y，Huang R Q，Densmore A L，etal.2009.Basic features and research progresses of Wenchuan Ms8.0earthquake.Journal of Sichuan University（in Chinese），41（3）：7-25.

Liu G X，Li J，Xu Z，etal.2010.Surface deformation associated with the 2008 Ms8.0Wenchuan earthquake from ALOS L-band SAR interferometry.International Journal of Applied Earth Observation and Geoinformation，12（6）：496-505，doi：10.1016／j.jag.2010.05.005.

Liu Z J，Zhang Z，Wen L，etal.2009.Co-seismic ruptures of the 12May 2008，Ms8.0Wenchuan earthquake，Sichuan：Eastwest crustal shortening on oblique，parallel thrusts along the eastern edge of Tibet.Earth Planet.Sci.Lett.，286（3-4）：355-370，doi：10.1016／j.epsl.2009.07.017.

Okada Y.1985.Surface deformation due to shear and tensile faults in a half-space.Bull.Seismol.Soc.Am，75（4）：1135-1154.

Shan X J，Qu C Y，Song X G.2009.Coseismic surface deformation caused by the Wenchuan Ms8.0earthquake from InSAR data analysis.Chinese Journal of Geophysics（in Chinese），52（2）：496-504.

Shen Z K，Sun J B，Zhang P Z，etal.2009.Slip maxima at fault junctions and rupturing of barriers during the 2008Wenchuan earthquake.Nature，10：718-724，doi：10.1038／ngeo636.

Tong X P，Sandwell D T，F(xiàn)ialko Y，etal.2010.Coseismic slip model of the 2008Wenchuan earthquake derived from joint inversion of interferometric synthetic aperture radar，GPS，and field data.Journal of Geophysical Research，115（B4）：1-19，doi：10.1029／2009JB006625.

Xu X W，Zhang P Z，Wen X Z，etal.2005.Features of active tectonics and recurrence behaviors of strong earthquakes in the western Sichuan province and its adjacent regions.Seismology and Geology（in Chinese），27（3）：446-461.

Zhang L，Wu J C，Chen Y L，etal.2007.Determination of coseismic displacement field of Taiwan Chi-Chi earthquake based on Doris.Journal of Geodesy and Geodynamics（in Chinese），27（3）：18-24.

Zhang P Z.2008.Co-seismic deformation field of 2008 Ms8.0 Wenchuan earthquake revealed by GPS data.Science in China Press（in Chinese），38（10）：1195-1206.

Zhang Y H，Gong W Y，Zhang J X，etal.2008.Measuring coseismic deformation of the sichuan earthquake by satellite differential INSAR.International Conference on Earth Observation Data Processing and Analysis，7285：4F-1：4F-8，doi：10.1117／12.816004.

Zhu A L，Xu X W，Diao G L，etal.2008.Relocation of the Ms8.0 Wenchuan earthquake sequence in part：preliminary seismotectonic analysis.Seismology and Geology（in Chinese），30（3）：759-767，doi：10.3969／j.issn.0253-4967.2008.03.014.

附中文參考文獻

朱艾斕，徐錫偉，刁桂苓等.2008.汶川Ms8.0地震部分余震重新定位及地震構造初步分析.地震地質(zhì)，30（3）：759-767，doi：10.3969／j.issn.0253-4967.2008.03.014.

張培震.2008.GPS測定的2008年汶川Ms8.0級地震的同震位移場.中國科學，38（10）：1195-1206.

張磊，伍吉倉，陳艷玲等.2007.基于Doris的臺灣集集地震同震位移場計算.大地測量與地球動力學，27（3）：18-24.

徐錫偉，張培震，聞學澤等.2005.川西及其鄰近地區(qū)活動構造基本特征與強震復發(fā)模型.地震地質(zhì)，27（3）：446-461.

單新建，屈春燕，宋小剛等.2009.汶川Ms8.0級地震InSAR同震形變場觀測與研究.地球物理學報，52（2）：496-504.

李勇，黃潤秋，Densmore A L等.2009.汶川8.0級地震的基本特征及其研究進展.四川大學學報.41（3）：7-25.

陳強，楊瑩輝，劉國祥等.2012（b）.基于邊界探測的InSAR最小二乘整周相位解纏方法.測繪學報，41（3）：441-448.

陳強，劉國祥，胡植慶等.2012（a）.GPS與PS-InSAR聯(lián)網(wǎng)監(jiān)測的臺灣屏東地區(qū)三維地表形變場.地球物理學報，55（10）：3248-3258，doi：10.6038／j.issn.0001-5733.2012.10.007.