大地水準(zhǔn)面高與insar形變測(cè)量不確定性的關(guān)系

大地水準(zhǔn)面高與insar形變測(cè)量不確定性的關(guān)系

1侵蝕、滑坡與海岸帶作為一種空間探測(cè)手段，聯(lián)合孔徑雷達(dá)干擾測(cè)量（iar）技術(shù)可以提供高精度、高空間分辨率、高區(qū)域地表變形數(shù)據(jù)。由于時(shí)間和幾何噪聲的顯著影響，以及氣候、道路、地表高度誤差和其他非結(jié)構(gòu)信號(hào)的顯著影響，傳統(tǒng)的iar方法對(duì)測(cè)量活動(dòng)斷層中的地震空間變形（如大區(qū)域、緩慢和小變形信號(hào)，主要限制）。自20世紀(jì)末以來(lái)，這兩種sar圖像時(shí)間序列分析算法（即永久散射體ps技術(shù)和短期差分收集bsa方法）的開發(fā)和開發(fā)不足，提高了測(cè)量形式變量的精度，并將其應(yīng)用于各種自然災(zāi)害、環(huán)境變化和人類活動(dòng)的表面變形過(guò)程，包括斷裂帶的地震損失。巖石和森林的積沙和顫抖；四川和冰蓋運(yùn)動(dòng)；造山運(yùn)動(dòng)和侵蝕、斜坡和海岸帶的變化；采礦、提取地下水、石油、天然氣等液體造成的礦區(qū)和城市地表沉降。重復(fù)軌道條帶模式(Strip-map/ImageMode)干涉圖受其幅寬限制(通常為100km×100km左右),難以捕獲到Mw>7大地震完整的同震形變場(chǎng),如2008年Mw7.9汶川地震、2011年Mw9.0日本Tohoku-Oki地震等,以及火山島弧、地殼內(nèi)部構(gòu)造應(yīng)變累積等引起的超大范圍形變場(chǎng).盡管通過(guò)拼接相鄰軌道數(shù)據(jù)方法可以進(jìn)行大范圍、連續(xù)的地表位移觀測(cè),但由于相鄰軌道數(shù)據(jù)獲取時(shí)間不同步,期間可能有多次強(qiáng)余震等噪聲干擾,因而相鄰軌道重合區(qū)經(jīng)平滑后不可避免會(huì)存在額外的誤差.與條帶模式相比,寬幅模式(WideSwathMode)SAR影像以降低空間分辨率為代價(jià)提供了更大的地面覆蓋,如EnvisatASAR寬幅模式數(shù)據(jù)的幅寬約為400km×400km,能夠?yàn)榉囱?Mw>7)同震破裂機(jī)制、震間滑動(dòng)速率與震后巖石圈流變機(jī)制等大范圍地殼運(yùn)動(dòng)與動(dòng)力學(xué)提供更完整的約束.目前通用的數(shù)字高程模型,如SRTMDEM和ASTERGDEM,都是基于EGM96大地水準(zhǔn)參考面,而常規(guī)InSAR干涉處理和高級(jí)InSAR時(shí)序分析等都需要在參考橢球面(如WGS84)下進(jìn)行.大地水準(zhǔn)面和參考橢球面之間的差異,即大地水準(zhǔn)面高(GeoidHeight)在全球不同區(qū)域分布都不一樣,其變化范圍約為[-106m,85m].絕大多數(shù)InSAR用戶都忽略了大地水準(zhǔn)面高,直接利用基于EGM96大地水準(zhǔn)參考面的SRTMDEM或ASTERDEM去除地形相位,大地水準(zhǔn)面高引起的相位誤差必然會(huì)對(duì)InSAR大范圍mm級(jí)地殼形變監(jiān)測(cè)應(yīng)用產(chǎn)生影響.另一方面,目前尚無(wú)公開發(fā)表的文獻(xiàn)定量分析大地水準(zhǔn)面高對(duì)InSAR形變監(jiān)測(cè)影響.因此,本文將以青藏高原北部阿爾金斷裂帶(AltynTaghFault)西段為例,討論InSAR數(shù)據(jù)處理中大地水準(zhǔn)面高對(duì)EnvisatASAR寬幅和條帶模式干涉處理結(jié)果的直接影響,并給出相應(yīng)的改正方法.2u3000egm98的或egm98的相位特征大地水準(zhǔn)面(Geoid)是與平均海水面最接近的重力等位面,而地球表面某一點(diǎn)的正高(Orthometricheight)是指該點(diǎn)沿鉛垂線至大地水準(zhǔn)面的距離.參考橢球(ReferenceEllipsoid)是數(shù)學(xué)上定義的地球表面,地球表面某一點(diǎn)的橢球高(Ellipsoidalheight)是指該點(diǎn)沿法線至參考橢球面的距離.橢球高Hellip與正高Hortho的關(guān)系如(1)式所示:Hellip=Hortho+N?(1)Ηellip=Ηortho+Ν?(1)其中,N代表橢球高與正高的差異,即大地水準(zhǔn)面高或大地水準(zhǔn)面差距(Geoidheight).WGS84參考橢球系統(tǒng)是為GPS全球定位系統(tǒng)使用而建立的坐標(biāo)系統(tǒng),同時(shí)也是空間大地測(cè)量(如InSAR等)與工程測(cè)量領(lǐng)域應(yīng)用較為廣泛的參考基準(zhǔn).理論上,若已知GPS測(cè)站的橢球高,基于現(xiàn)存的大地水準(zhǔn)面模型(如EGM96)就可以確定相應(yīng)的正高.WGS84橢球高與EGM96正高系統(tǒng)之間的轉(zhuǎn)換關(guān)系如(2)式所示:HWGS84=HEGM96+N.(2)ΗWGS84=ΗEGΜ96+Ν.(2)本文略去推導(dǎo)過(guò)程,直接給出重復(fù)軌道干涉測(cè)量相位的一般表達(dá)式:?=?ref+?topo+?defo+?noise=?4πλ(B//+B⊥RsinθHWGS84?Δρ)+?noise,(3)?=?ref+?topo+?defo+?noise=-4πλ(B//+B⊥RsinθΗWGS84-Δρ)+?noise,(3)其中θ為入射角(incidenceangle),R為衛(wèi)星至地面點(diǎn)的距離,?ref為參考相位,即地球曲率相位貢獻(xiàn),?topo為地形相位,?defo為形變相位,?noise為大氣、軌道與熱噪聲等相位貢獻(xiàn).若給定了外部DEM數(shù)據(jù),就可以模擬SAR影像得到地形相位?topo.理想情況下忽略噪聲相位?noise,根據(jù)式(3)可以由去除參考相位?ref之后的展平相位?flat和地形相位?topo求得形變量Δρ,如下所示:Δρ=λ4π?defo=λ4π(?flat??topo)=λ4π?flat+B⊥RsinθHWGS84.(4)Δρ=λ4π?defo=λ4π(?flat-?topo)=λ4π?flat+B⊥RsinθΗWGS84.(4)根據(jù)誤差傳播定律,由(4)式可知外部DEM中的地形誤差與形變測(cè)量不確定性成正比關(guān)系,即σΔρ=B⊥RsinθσDEM.(5)σΔρ=B⊥RsinθσDEΜ.(5)如果使用的外部DEM基于EGM96正高系統(tǒng),由(2)式和(3)式可知由大地水面高引起地形相位誤差為?N,即EGM96大地水準(zhǔn)面高在WGS84參考橢球系下的相位貢獻(xiàn):?N=?topo??EGM96=(?4πλB⊥RsinθHWGS84)?(?4πλB⊥RsinθHEGM96)=?4πλB⊥RsinθN.(6)?Ν=?topo-?EGΜ96=(-4πλB⊥RsinθΗWGS84)-(-4πλB⊥RsinθΗEGΜ96)=-4πλB⊥RsinθΝ.(6)3實(shí)驗(yàn)與討論3.1基于gps的dem誤差測(cè)量阿爾金斷裂帶(AltynTaghFault)位于青藏高原北部,北臨塔里木盆地,南臨青藏高原主板塊,是該地區(qū)最顯著(～2000km)的北東東向左旋走滑斷裂帶,由一系列呈羽狀、雁列分布的斷裂組合而成,是研究?jī)?nèi)陸地震周期性質(zhì)與大陸構(gòu)造變形機(jī)制的天然實(shí)驗(yàn)室.在其最西端,阿爾金斷裂帶與康西瓦左旋走滑斷裂帶(KarakaxFault)和西昆侖逆沖帶(WesternKunlunThrustBelt)相連;在其西南端,阿爾金斷裂帶與相鄰的龍木—郭扎錯(cuò)左旋走滑斷裂帶(Longmu-GozhaCoFault)相連.本文的研究區(qū)域主要覆蓋阿爾金斷裂帶西段和龍木—郭扎錯(cuò)斷裂帶.如圖1a所示,紅色方框代表本文所使用的EnvisatASAR寬幅影像覆蓋范圍,藍(lán)色矩形代表?xiàng)l帶模式影像的覆蓋范圍.圖1b給出了由塔里木盆地穿越阿爾金斷裂帶至青藏高原做地形剖面線的結(jié)果,高程從北到南由1000m左右陡升至6000m左右.大地水準(zhǔn)面和重力異常的劇烈橫向變化均與陡峭的地形梯度有關(guān).目前SRTMDEM和ASTERGDEM是InSAR干涉處理中常用的數(shù)字高程模型.本文采用2009—2011年中國(guó)陸態(tài)環(huán)境網(wǎng)絡(luò)工程GPS區(qū)域站點(diǎn)坐標(biāo)驗(yàn)證了這兩種DEM在該地區(qū)的高程質(zhì)量(圖2a).中國(guó)陸態(tài)環(huán)境網(wǎng)絡(luò)工程GPS數(shù)據(jù)的三維坐標(biāo)精度優(yōu)于厘米級(jí),并已由WGS84橢球坐標(biāo)系轉(zhuǎn)換至EGM96正高系統(tǒng).如圖2b和2c所示,ASTERGDEM-GPS的誤差均值為0.1m,而SRTMDEM-GPS的誤差均值為-5.6m;但SRTM-GPS的高程精度(8.8m)略優(yōu)于ASTERGDEM-GPS(9.6m),這與已發(fā)表的結(jié)果基本一致.本文采用SRTM作為InSAR干涉處理的外部高程數(shù)據(jù)源,根據(jù)式(5)對(duì)DEM誤差與寬幅、條帶模式InSAR形變測(cè)量不確定性關(guān)系進(jìn)行模擬.利用JPL/Caltech開發(fā)的ROI_PAC開源軟件和北京大學(xué)遙感所開發(fā)的寬幅InSAR程序處理了EnvisatASAR條帶模式(Track205共44景:2003-04—2010-08)和寬幅模式(Track434共27景:2007-02—2010-09)原始數(shù)據(jù),分別生成946和351幅干涉影像.圖3a和3b分別給出了寬幅模式干涉對(duì)(Track434)和條帶模式干涉對(duì)(Track205)對(duì)應(yīng)的空間垂直基線分布的直方圖;由于本文處理的寬幅影像均為2007年1月23日Envisat衛(wèi)星軌道優(yōu)化調(diào)整后的數(shù)據(jù),寬幅干涉影像的垂直基線都在900m以內(nèi)(圖3a),而條帶模式包含衛(wèi)星軌道調(diào)整前的數(shù)據(jù),其干涉影像的垂直基線最長(zhǎng)達(dá)2000m(圖3d).由于不同波長(zhǎng)對(duì)應(yīng)不同的臨界基線值Bcritical,超過(guò)該值可能導(dǎo)致相干性的完全損失.對(duì)于L波段的數(shù)據(jù),其臨界基線值可達(dá)6.5～13.1km;對(duì)于EnvisatASAR(C波段)數(shù)據(jù),其寬幅模式和條帶模式干涉測(cè)量的臨界基線值約為450m和1100m.由圖3a和圖3d可知,均有90%的寬幅和條帶模式干涉對(duì)垂直基線小于上述臨界基線值,因此在本文接下來(lái)模擬分析中可以忽略臨界基線值的約束,即考慮所有干涉對(duì)的基線值.由圖2c可知該地區(qū)的SRTMDEM的高程誤差為8.8m,利用公式(5)可以計(jì)算得到所有InSAR干涉影像中DEM所引起的誤差統(tǒng)計(jì)直方圖(圖3b和3e).對(duì)于InSAR時(shí)間序列分析而言,我們通常希望得到優(yōu)于1.2mm/yr的精度,因此利用時(shí)間跨度近4年的寬幅數(shù)據(jù)測(cè)定的形變量的誤差應(yīng)低于5mm,時(shí)間跨度為7年的條帶模式數(shù)據(jù)測(cè)定的形變量的誤差應(yīng)不高于8mm.由圖3b和3e可知,約有49%寬幅干涉對(duì)的地形所引起的誤差超過(guò)5mm,約有64%條帶模式干涉對(duì)的地形所引起的誤差超過(guò)8mm.由圖3c和3f可知空間基線的大小決定了InSAR干涉圖對(duì)地形誤差的敏感性.當(dāng)DEM誤差為8.8m,空間基線為30～100m,相應(yīng)地會(huì)引入0.7～3mm的形變測(cè)量不確定性;當(dāng)空間基線為1100m(最大相干基線)時(shí)將會(huì)引入超過(guò)25mm的誤差.3.2巖石圈結(jié)構(gòu)的特殊及地殼水準(zhǔn)面特征EMG96大地水準(zhǔn)面高(N)代表了由EGM96重力等位面確定的大地水準(zhǔn)面與WGS84參考橢球的差異,其在全球的分布如圖4a所示,青藏高原、印度、西非、南美洲南部、北美洲北部等地區(qū)存在明顯的負(fù)大地水準(zhǔn)面高,印度洋中部具有最大的大地水準(zhǔn)面高負(fù)值(約-106m),而其他地區(qū)如印尼、菲律賓、澳大利亞、歐洲、非洲和中美洲的大部分地區(qū)存在正的大地水準(zhǔn)面高.青藏高原北邊緣阿爾金斷裂帶地區(qū)具有顯著的高程、重力與大地水準(zhǔn)面異常變化,與其特殊的巖石圈結(jié)構(gòu)密切相關(guān).位于阿爾金斷裂帶西段以南,與喀喇昆侖斷裂帶包圍的青藏高原主板塊的EGM96大地水準(zhǔn)面高在-20～-30m左右;而位于阿爾金斷裂帶以北的塔里木盆地和喀喇昆侖斷裂帶以南的區(qū)域?yàn)?50～-60m左右(圖4b和4c).張赤軍等(2011)認(rèn)為塔里木盆地和剛蓋盆地下部深遠(yuǎn)處異常物質(zhì)效應(yīng)是導(dǎo)致青藏高原EGM96大地水準(zhǔn)面高存在較大負(fù)值分布的主要原因,而塔里木盆地及以南莫霍面的突變直接加劇了昆侖山大地水準(zhǔn)面的梯度變化.其他學(xué)者認(rèn)為對(duì)地幔結(jié)構(gòu)具有強(qiáng)約束的動(dòng)力地形(DynamicTopography)導(dǎo)致的質(zhì)量異常對(duì)該地區(qū)大地水準(zhǔn)面有較大影響.如圖4b與4c所示,本文所用寬幅數(shù)據(jù)434軌道所覆蓋區(qū)域的大地水準(zhǔn)面高在南北方向存在顯著的梯度變化,變化區(qū)間在[-23m,-53m],即變化范圍達(dá)到了30m.同樣,對(duì)于條帶模式數(shù)據(jù)205軌道所覆蓋區(qū)域變化量也達(dá)到了20m.若將大地水準(zhǔn)面高變化量考慮為大地水準(zhǔn)面高可能引進(jìn)的誤差,由圖3c和3f可知,當(dāng)垂直基線為100m時(shí)將對(duì)條帶模式和寬幅模式InSAR形變測(cè)量分別產(chǎn)生7.6mm和8.8mm的影響.3.3dem的建立本文僅分別選取寬幅和條帶模式時(shí)間基線相對(duì)較短的3個(gè)干涉對(duì),作為實(shí)驗(yàn)數(shù)據(jù)來(lái)考察該地區(qū)大地水準(zhǔn)面高的變化對(duì)InSAR形變測(cè)量的影響.除受空間去相關(guān)和時(shí)間去相關(guān)的影響之外,寬幅干涉測(cè)量還需要主輔影像對(duì)之間的方位向掃描模式保持同步,其同步性比例(AzimuthScanningPatternSynchronization/BurstOverlap)的高低直接決定了寬幅干涉圖相干性的好壞.如表1和表2所示,3個(gè)寬幅干涉對(duì)和3個(gè)條帶模式干涉對(duì)基本處于相同季節(jié),空間基線≤125m,時(shí)間基線≤70天,寬幅數(shù)據(jù)干涉對(duì)的信號(hào)片段(burst)同步性比例均高于60%.采用歐洲空間局DORIS精密軌道(VOR)用于空間基線計(jì)算,并使用分辨率為3″的SRTMversion4.1數(shù)據(jù)去除地形相位,采用功率譜濾波以提高信噪比,并使用SNAPHU方法對(duì)濾波之后干涉圖進(jìn)行解纏,如圖5(a1,b1,c1)和圖6(a1,b1,c1)所示.本文采用以下方法對(duì)InSAR數(shù)據(jù)處理中的DEM進(jìn)行大地水準(zhǔn)面高改正,即將其高程系統(tǒng)由EGM96模型轉(zhuǎn)換為WGS84參考橢球系:(1)首先,需要將研究區(qū)域的大地水準(zhǔn)面高與DEM的分辨率保持一致.本文采用IDL程序雙線性內(nèi)插方法將EGM96大地水準(zhǔn)面高由15′×15′插值到3″×3″的分辨率.(2)根據(jù)公式(2)將SRTM的EGM96高程改正至WGS84參考系.考慮到海洋上的大地水準(zhǔn)面高與通常默認(rèn)的零海平面高可能存在較大差異,如印度洋的中部和太平洋的馬里亞納海溝地區(qū).為統(tǒng)一起見,本文將改正后的DEM仍然假定其海平面區(qū)域?yàn)榱愀叱?注意SAR影像在海洋上的干涉測(cè)量應(yīng)用不在本文的研究范圍.3.4不確定區(qū)域內(nèi)相位差的消除利用改正后的DEM生成新的一組差分干涉測(cè)量結(jié)果,如圖5(a2,b2,c2)和圖6(a2,b2,c2)所示.圖5(a3,b3,c3)顯示了利用大地水準(zhǔn)面高改正前后的SRTMDEM生成的寬幅模式干涉影像對(duì)的相位差值,即代表大地水準(zhǔn)面高引入的誤差;由圖5(a3,b3,c3)可見,隨干涉對(duì)空間垂直基線的增加,干涉圖中殘余相位(?N)的條紋數(shù)亦隨之增加.根據(jù)公式(6)可知,該殘余相位的變化量與垂直基線的變化量成正比.對(duì)3個(gè)條帶模式數(shù)據(jù)干涉對(duì)進(jìn)行與寬幅模式數(shù)據(jù)類似的檢驗(yàn),由圖6(c3)可知存在較為明顯的相位差,與條帶模式數(shù)據(jù)所在區(qū)域的大地水準(zhǔn)面高的趨勢(shì)相同.圖6(b3,c3)要比圖6(a3)對(duì)應(yīng)的干涉對(duì)垂直基線大,由此引入了更大的大地水準(zhǔn)面高相位?N.基于空間大地測(cè)量觀測(cè)手段(含GPS和InSAR)的研究顯示阿爾金斷裂帶西段的平均滑動(dòng)速率<10mm/yr,屬于大范圍、微小形變的范疇.在InSAR干涉影像的解析過(guò)程中,不少用戶通常利用最佳擬合平面或者曲面去削弱衛(wèi)星軌道誤差的影響.然而,對(duì)于長(zhǎng)波段地殼形變信號(hào)而言,去除最佳擬合曲面以改正軌道誤差的方法很大程度上會(huì)同時(shí)去除長(zhǎng)波段的形變信號(hào).另一方面,本文也嘗試對(duì)DEM改正前后得到的干涉圖相位差去除最佳擬合平面,得到了去平面后的相位差的殘差,如圖5(a4,b4,c4)和圖6(a4,b4,c4)所示.(1)顯而易見,當(dāng)空間垂直基線較長(zhǎng)(如:B>70m)時(shí),去除最佳擬合平面的方法也難以消除大地水準(zhǔn)面高引進(jìn)的誤差(圖5(b4,c4)和圖6(c4)).由于大地水準(zhǔn)面高梯度變化趨勢(shì)與該區(qū)域的地形存在較大相關(guān)性,而與軌道條紋模式不相關(guān),當(dāng)使用InSAR測(cè)量長(zhǎng)波段地殼形變時(shí),通過(guò)去除軌道誤差的方法并不適合、也不可能完全改正大地水準(zhǔn)面高的變化的影響.(2)由EGM96大地水準(zhǔn)面高的全球分布特征(圖4a)可以看出,大地水準(zhǔn)面高的影響并非僅與研究范圍大小有關(guān),同時(shí)也與研究區(qū)域位置相關(guān).對(duì)于大地水準(zhǔn)面高變化比較緩慢的區(qū)域,若采用短基線(如B<70m)干涉圖、去除二次項(xiàng)軌道平面的方法一定程度上可以削弱因大地水準(zhǔn)面高引入的誤差.但是,短基線的條件有些情況很難滿足.綜上所述,當(dāng)利用InSAR方法測(cè)定全球大陸地區(qū)的大范圍(如研究范圍大于100km×100km)、毫米級(jí)微小形變時(shí),使用基于WGS84參考橢球體的DEM應(yīng)該是首選.4對(duì)地表變形的影響由于InSAR數(shù)據(jù)處理軟件(如ROI_PAC等)采用WGS84橢球高坐標(biāo)系統(tǒng),而常規(guī)使用的DEM(如SRTM或ASTERGDEM)均為EGM96正高系統(tǒng),本文評(píng)估了不同高程系統(tǒng)對(duì)EnvisatASAR寬幅和條帶模式干涉測(cè)量在阿爾金斷裂帶西段地殼形變測(cè)量精度的影響,得到以下幾點(diǎn)結(jié)論:(1)由理論模擬分析顯示,對(duì)于100m的垂直基線,8.8m的SRTMDEM測(cè)量誤差,若存在2