DEM輔助偏移量跟蹤技術(shù)的山地冰川運(yùn)動(dòng)監(jiān)測(cè)研究

王 群， 范景輝， 周 偉， 袁蔚林， 童立強(qiáng)， 郭兆成

(1.中國(guó)地質(zhì)大學(xué)(北京)土地科學(xué)技術(shù)學(xué)院，北京 100083；2.中國(guó)國(guó)土資源航空物探遙感中心，北京 100083)

0 引言

山地冰川是冰川的重要組成部分，是反映全球氣候變化的敏感的記錄器和指示器[1]。冰川運(yùn)動(dòng)是冰川的重要特征之一，能引起許多重大的自然災(zāi)害[2]。其中，冰川湖潰決及其引起的洪水災(zāi)害常給受影響區(qū)帶來巨大的危害[3-4]。因此，監(jiān)測(cè)山地冰川的運(yùn)動(dòng)狀況對(duì)全球氣候變化研究和冰川湖潰決危險(xiǎn)性評(píng)估十分重要[5-6]。

合成孔徑雷達(dá)(synthetic aperture Radar, SAR)可以全天時(shí)、全天候獲取數(shù)據(jù)，能夠長(zhǎng)時(shí)間、大范圍地監(jiān)測(cè)山地冰川區(qū)地表環(huán)境[7]，相較受到時(shí)空限制的傳統(tǒng)實(shí)地監(jiān)測(cè)[8]和受到多云雨天氣限制的光學(xué)遙感[9]更具優(yōu)勢(shì)，是冰川研究中的重要工具[10]。過去30多a以來，已有學(xué)者使用基于SAR影像的差分干涉測(cè)量技術(shù)(differential interferometry synthetic aperture Radar，DInSAR)、偏移量跟蹤技術(shù)以及多孔徑雷達(dá)干涉測(cè)量技術(shù)(multiple aperture InSAR，MAI)成功地監(jiān)測(cè)了冰川表面運(yùn)動(dòng)情況[10-12]。但是，DInSAR技術(shù)和MAI技術(shù)主要利用SAR相位信息，對(duì)相干性要求較高，而冰川流速較快或SAR影像空間基線較大、時(shí)間基線較長(zhǎng)時(shí)往往會(huì)使冰川區(qū)域失相干，限制了DInSAR和MAI技術(shù)的應(yīng)用[9,13]。基于SAR影像強(qiáng)度信息的偏移量跟蹤技術(shù)，在SAR影像之間不能保持良好相干性的情況下，仍可以直接獲得冰川在方位向和距離向的位移[14]，在監(jiān)測(cè)冰川運(yùn)動(dòng)方面比其他2種技術(shù)適用范圍更大。

山地冰川大多分布于高海拔、地形陡峭的區(qū)域，傳統(tǒng)的偏移量跟蹤技術(shù)為了減小地形的影響，常采用空間基線較小的SAR影像數(shù)據(jù)來進(jìn)行冰川運(yùn)動(dòng)速度的測(cè)量[15]，這大大限制了數(shù)據(jù)選取范圍，不利于滿足冰川研究需要。Yan等[16]針對(duì)地形崎嶇地區(qū)、大空間基線SAR影像條件下地形引起的偏移量誤差，提出一種基于外部數(shù)字高程模型(digital elevation model，DEM)的地形影響補(bǔ)償算法，通過計(jì)算出一定基線條件下地形起伏引起的方位向和距離向偏移量，并從總的偏移量中予以去除，從而獲得更為精確的結(jié)果?；诖怂惴ㄊ褂肔波段PALSAR數(shù)據(jù)在西昆侖山東部冰川區(qū)[16]、慕士塔格峰冰川區(qū)[17]以及塔吉克斯坦南部冰川區(qū)[18]獲取了精確的冰川表面運(yùn)動(dòng)信息，證明了算法的可靠性。本文DEM輔助的偏移量跟蹤技術(shù)利用DEM和SAR影像軌道信息建立起主從影像的配準(zhǔn)查找表，將研究區(qū)地形引起的偏移量計(jì)算到查找表內(nèi)，在空間基線較大時(shí)能更為精確地配準(zhǔn)2幅影像，不僅在監(jiān)測(cè)山地冰川表面運(yùn)動(dòng)時(shí)有效打破了基線長(zhǎng)度的限制，也能應(yīng)用于同一傳感器不同波束模式的SAR影像，進(jìn)一步擴(kuò)大了監(jiān)測(cè)冰川運(yùn)動(dòng)所使用的SAR數(shù)據(jù)范圍。

本文使用2對(duì)不同空間基線的TerraSAR-X數(shù)據(jù)和一對(duì)不同入射角的COSMO-SkyMed數(shù)據(jù)，以卓莫拉日山系東段冰川分布區(qū)為研究區(qū)，對(duì)TerraSAR-X應(yīng)用傳統(tǒng)的以及DEM輔助的偏移量跟蹤技術(shù)獲取了研究區(qū)方位向和距離向的偏移量分布情況，通過對(duì)比2種結(jié)果在非冰川區(qū)的偏移量分布情況以及均方根誤差，顯示了本文方法的優(yōu)勢(shì)； 對(duì)COSMO-SkyMed進(jìn)行數(shù)據(jù)處理獲取了冰川表面運(yùn)動(dòng)信息，表明了本文方法可以進(jìn)一步擴(kuò)大可用的SAR影像來源。

1 技術(shù)方法

1.1 地形引起偏移量分析

偏移量跟蹤技術(shù)利用歸一化互相關(guān)算法來獲取SAR影像方位向和距離向偏移量[8]，其原理是在主影像中選取一塊特定的窗口，并將它和從影像中的搜索窗口進(jìn)行匹配，并計(jì)算每一個(gè)位置的互相關(guān)系數(shù)，當(dāng)互相關(guān)系數(shù)值最大時(shí)，匹配完成，從而獲得匹配位置的總偏移量。

在地形起伏不大、基線較小，且不考慮大氣擾動(dòng)影響的情況下，總偏移量包括冰川位移引起的偏移量和軌道偏移量[19]。傳統(tǒng)的偏移量跟蹤技術(shù)在圖像配準(zhǔn)階段基于圖像亞像元互相關(guān)算法計(jì)算2幅影像的同名點(diǎn)，尋找影像的幾何變換關(guān)系，一般使用最小二乘法計(jì)算主從影像的多項(xiàng)式變換模型，實(shí)現(xiàn)2幅影像的配準(zhǔn)[20]。但在地形變化劇烈以及大空間基線條件下，地形引起的SAR影像之間的局部配準(zhǔn)偏移量不可忽視[21]。傳統(tǒng)的偏移量跟蹤技術(shù)忽略了地形引起的偏移量，很難對(duì)受此影響的SAR影像進(jìn)行精確配準(zhǔn)，獲得準(zhǔn)確的運(yùn)動(dòng)監(jiān)測(cè)結(jié)果。

地形在方位向和距離向引起的偏移量為[22]

(1)

(2)

式中：Atopo和Rtopo分別為地形引起的方位向偏移量和距離向偏移量；α為軌道交角；θ為入射角；Δh為高程差；a和r分別為方位向和距離向的分辨率；B⊥為垂直基線長(zhǎng)度；R為斜距長(zhǎng)。

根據(jù)以上公式可以看出，方位向偏移量與地形起伏、軌道交角以及入射角等有關(guān)； 在距離向上地形起伏導(dǎo)致的偏移量大小與垂直基線成正比。為了更直觀地了解2分量與地形起伏之間的關(guān)系，參照實(shí)際情況設(shè)定TerraSAR-X數(shù)據(jù)θ為26°，r為0.9 m，a為2.0 m，R為560 km。則由式(1)可知，若主從影像對(duì)應(yīng)的軌道交角α為0.025°，與平地條件相比，280 m的地形起伏將為同名地物在主從影像方位向上帶來1/8個(gè)像元的額外偏移； 當(dāng)然，若主從影像對(duì)應(yīng)的軌道平行，地形起伏將不會(huì)引入相應(yīng)偏移。由式(2)可知，100 m長(zhǎng)度垂直基線情況下，當(dāng)研究區(qū)地形起伏達(dá)到300 m時(shí)，將為同名地物在主從影像距離向上帶來1/8個(gè)像元的偏移； 1 000 m長(zhǎng)度垂直基線情況下，研究區(qū)地形起伏達(dá)到30 m時(shí)，將為同名地物在主從影像距離向上帶來1/8個(gè)像元的偏移。因此，大空間基線條件下，地形起伏引起的距離向偏移量已經(jīng)不可忽略，傳統(tǒng)偏移量跟蹤技術(shù)忽略了地形因素，大空間基線時(shí)無法對(duì)SAR影像精確配準(zhǔn)，運(yùn)動(dòng)監(jiān)測(cè)結(jié)果的精度會(huì)受到較大影響。

1.2 DEM輔助的偏移量跟蹤技術(shù)

DEM輔助的偏移量跟蹤技術(shù)，核心在于利用研究區(qū)的DEM和衛(wèi)星軌道信息建立2幅SAR影像之間的映射關(guān)系，即配準(zhǔn)查找表。因?yàn)樵摷夹g(shù)將地形引起的偏移量也計(jì)算在查找表內(nèi)，所以可以獲得更為精確的配準(zhǔn)[23]。在地形變化劇烈及大空間基線條件下DEM輔助的偏移量跟蹤處理流程如圖1所示。主要分為3個(gè)步驟：

1)初始配準(zhǔn)查找表的生成。地理編碼獲得基于主影像成像幾何的SAR坐標(biāo)系DEM及兩者之間的映射關(guān)系； 然后建立起SAR坐標(biāo)系DEM與從影像的映射關(guān)系； 結(jié)合2種映射關(guān)系，生成了記錄主從影像偏移量的初始配準(zhǔn)查找表，該查找表已將系統(tǒng)偏移與地形引起的偏移量同時(shí)計(jì)算在內(nèi)。

2)初始配準(zhǔn)查找表的精化。利用重采樣獲得具有相同成像幾何的主、從影像，進(jìn)一步使用圖像互相關(guān)方法獲取兩者的配準(zhǔn)多項(xiàng)式，據(jù)此精化配準(zhǔn)查找表，降低可能存在的由軌道定位和DEM數(shù)據(jù)不精確帶來的誤差。

3)位移場(chǎng)的計(jì)算。依據(jù)精化配準(zhǔn)查找表重采樣從影像，按照傳統(tǒng)偏移量跟蹤技術(shù)流程計(jì)算出冰川方位向和距離向位移場(chǎng)。

冀東油田勘探開發(fā)區(qū)塊分為南堡灘海和南堡陸地，鉆遇中淺層地層主要為平原組、明化鎮(zhèn)組、館陶組和東營(yíng)組。其中，南堡陸地部分區(qū)塊東營(yíng)組缺失，沙河街地層溫度相對(duì)較低，約110℃左右;上部平原組、明化鎮(zhèn)組和館陶組多以砂泥巖為主，黏土含量高，鉆井過程中巖屑水化分散嚴(yán)重，快鉆時(shí)要求鉆井液體系有較強(qiáng)的包被抑制能力，館陶組底部普遍發(fā)育玄武巖和底礫巖，要求鉆井液體系有較好的抑制能力和封堵能力。鉆遇儲(chǔ)層要求鉆井液體系具有良好的降濾失性能和高中孔滲儲(chǔ)層保護(hù)性能。

此方法能夠獲得更精確的偏移量結(jié)果，對(duì)同一傳感器不同波束模式的SAR影像也適用。

圖1 外部DEM輔助偏移量跟蹤技術(shù)處理流程圖Fig.1 Flow chart of the external DEM-assisted offset tracking technique

2 研究區(qū)概況與數(shù)據(jù)源

2.1 研究區(qū)概況

研究區(qū)位于西藏自治區(qū)康馬縣和浪卡子縣之間喜馬拉雅山脈東部、卓莫拉日山系東段，N28°07′～28°19′，E89°53′～90°22′之間，地處中國(guó)和不丹交界處，地形起伏較大，高程差超過2 000 m，平均海拔在5 000 m以上，包括大小近10條冰川。同時(shí)，研究區(qū)內(nèi)有桑旺錯(cuò)和黃湖2個(gè)較大的冰川湖，其中桑旺錯(cuò)于1954年發(fā)生潰決，受災(zāi)人數(shù)達(dá)2萬余人，死亡人數(shù)約400人，造成巨大損失[24]，現(xiàn)今仍有重大的成災(zāi)風(fēng)險(xiǎn)[25]。因此，監(jiān)測(cè)研究區(qū)的冰川運(yùn)動(dòng)情況對(duì)冰川湖潰決危險(xiǎn)評(píng)估和當(dāng)?shù)厣a(chǎn)生活安全有重要意義。

2.2 數(shù)據(jù)來源與處理

TerraSAR-X 雷達(dá)衛(wèi)星具有多極化、多入射角和精確的姿態(tài)和軌道控制能力,采用波長(zhǎng)約3.2 cm的X波段SAR,可以全天時(shí)、全天候地對(duì)地觀測(cè)[26]。COSMO-SkyMed對(duì)地觀測(cè)系統(tǒng)是由4顆搭載了多個(gè)成像模式的X波段SAR的太陽同步軌道衛(wèi)星組成，在軌道上相鄰的不同衛(wèi)星對(duì)同一地點(diǎn)的重訪時(shí)間分別為1 d,3 d,4 d和8 d[27]。

選取2對(duì)不同長(zhǎng)度空間基線、時(shí)間間隔為11 d的TerraSAR-X數(shù)據(jù)，主要成像參數(shù)為升軌右視、條帶掃描模式，空間分辨率為3 m，其數(shù)據(jù)參數(shù)如表1所示。選取一對(duì)不同波束模式的COSMO-SkyMed影像，為降軌右視獲取，空間分辨率為3 m，波束模式分別為HI01和HI02，2景影像在入射角和覆蓋范圍上均有差異，數(shù)據(jù)參數(shù)如表2所示。DEM為30 m空間分辨率的SRTM，絕對(duì)高程精度約為15 m[28]。研究區(qū)范圍及數(shù)據(jù)覆蓋范圍如圖2所示。

表1 TerraSAR-X數(shù)據(jù)參數(shù)Tab.1 Parameters of TerraSAR-X data

表2 COSMO-SkyMed數(shù)據(jù)參數(shù)Tab.2 Parameters of COSMO-SkyMed data

圖2 研究區(qū)范圍及TerraSAR-X數(shù)據(jù)、不同入 射角COSMO-SkyMed數(shù)據(jù)覆蓋范圍Fig.2 Study area and TerraSAR-X covering map and the covering map of different incidence COSMO-SkyMed data

為了提高影像處理效率，將TerraSAR-X和COSMO-SkyMed影像范圍適當(dāng)裁剪。經(jīng)過多次試驗(yàn)發(fā)現(xiàn)，偏移量跟蹤搜索窗口過大時(shí)，不利于保存細(xì)節(jié)信息，而搜索窗口過小則會(huì)使冰川的位移情況過于細(xì)碎，不利于分析，故本文將偏移量跟蹤的搜索窗口設(shè)為100像元×100像元，既可以較好地保留細(xì)節(jié)，又使得位移信息不太零散，同時(shí)也能防止遺漏偏移量較大的點(diǎn)。

3 結(jié)果分析

3.1 偏移量跟蹤技術(shù)結(jié)果

對(duì)TerraSAR-X應(yīng)用傳統(tǒng)及DEM輔助的偏移量跟蹤技術(shù)獲取的非冰川區(qū)監(jiān)測(cè)結(jié)果如圖3和圖4所示。圖例中不同顏色表示偏移量值的變化，0值代表無偏移量，每一個(gè)色周代表3 m的偏移量，當(dāng)研究區(qū)某處偏移量絕對(duì)值大于3 m時(shí)，偏移量結(jié)果圖上將顯示出色棒中顏色重復(fù)出現(xiàn)構(gòu)成的顏色條紋。黑色區(qū)域?yàn)橄嚓P(guān)性閾值低于設(shè)定閾值的區(qū)域及SAR成像時(shí)的疊掩陰影區(qū)，即無結(jié)果區(qū)域。

(a) 傳統(tǒng)偏移量跟蹤技術(shù)方位向偏移量結(jié)果(b) DEM輔助偏移量跟蹤技術(shù)方位向偏移量結(jié)果

(c) 傳統(tǒng)偏移量跟蹤技術(shù)距離向偏移量結(jié)果(d) DEM輔助偏移量跟蹤技術(shù)距離向偏移量結(jié)果

圖320160625—20160706像對(duì)提取結(jié)果

Fig.3Offsettrackingresultsfromimagepair20160625-20160706

(a) 傳統(tǒng)偏移量跟蹤技術(shù)方位向偏移量結(jié)果 (b) DEM輔助偏移量跟蹤技術(shù)方位向偏移量結(jié)果

(c) 傳統(tǒng)偏移量跟蹤技術(shù)距離向偏移量結(jié)果 (d) DEM輔助偏移量跟蹤技術(shù)距離向偏移量結(jié)果

圖420160808—20160819像對(duì)提取結(jié)果

Fig.4Offsettrackingresultsfromimagepair20160808-20160819

從圖上可以直觀地看到DEM輔助的偏移量跟蹤結(jié)果在非冰川區(qū)接近0值的顏色更加干凈統(tǒng)一，表明這種方法所獲得的冰川偏移量結(jié)果更加可靠。從傳統(tǒng)方法在非冰川區(qū)方位向的偏移量監(jiān)測(cè)結(jié)果來看，隨著垂直基線的增加，非冰川區(qū)小部分區(qū)域誤差雖略有增加，但整體仍趨于0值，說明方位向上地形起伏造成的偏移量誤差受空間基線變化影響較??； 與之對(duì)比，DEM輔助的方法對(duì)監(jiān)測(cè)結(jié)果有所改善。從傳統(tǒng)方法在非冰川區(qū)距離向的偏移量監(jiān)測(cè)結(jié)果來看，隨著垂直基線的增大，非冰川區(qū)距離向的偏移量誤差更加明顯，甚至?xí)^一個(gè)像元大小，說明距離向上地形起伏造成的偏移量誤差隨著空間基線的增大而變大。與之對(duì)比可見，在空間基線較大時(shí)，DEM輔助的方法對(duì)距離向上非冰川區(qū)監(jiān)測(cè)結(jié)果的改善效果非常明顯。

對(duì)不同入射角、覆蓋范圍的COSMO-SkyMed應(yīng)用DEM輔助的偏移量跟蹤技術(shù)的監(jiān)測(cè)結(jié)果如圖5所示。由圖5可知，監(jiān)測(cè)結(jié)果在非冰川區(qū)接近0值，表明此方法同樣適用于不同入射角、覆蓋范圍的SAR數(shù)據(jù)。

(a) DEM輔助偏移量跟蹤技術(shù)方位向偏移量結(jié)果 (b) DEM輔助偏移量跟蹤技術(shù)距離向偏移量結(jié)果

圖520160630-20160731像對(duì)提取結(jié)果

Fig.5Offsettrackingresultsfromimagepair20160630-20160731

3.2 誤差分析

為了評(píng)價(jià)山地冰川流速估算的精度，對(duì)TerraSAR-X和COSMO-SkyMed結(jié)果非冰川區(qū)分別選取一塊區(qū)域，計(jì)算區(qū)域內(nèi)位移殘差的均方根誤差(root mean square error,RMSE)來衡量結(jié)果的精度。對(duì)TerraSAR-X數(shù)據(jù)處理結(jié)果選取的區(qū)域范圍如圖3(a)紅色矩形所示，RMSE如表3所示； 對(duì)COSMO-SkyMed數(shù)據(jù)處理結(jié)果選取的區(qū)域范圍如圖5(a)紅色矩形所示，RMSE如表4所示。

表3TerraSAR-X數(shù)據(jù)傳統(tǒng)和DEM輔助的

偏移量跟蹤技術(shù)非冰川區(qū)RMSE

Tab.3RMSEinice-freeregionofTraditionalandexternalDEM-assistedoffsettrackingtechniquesusingTerraSAR-Xdata

(m)

數(shù)據(jù)對(duì)傳統(tǒng)偏移量跟蹤技術(shù)DEM輔助的偏移量跟蹤技術(shù)方位向偏移量RMSE距離向偏移量RMSE方位向偏移量RMSE距離向偏移量RMSE20160625/201607060.0470.0740.0410.03620160808/201608190.1091.9880.0930.089

表4不同入射角COSMO-SkyMed數(shù)據(jù)

偏移量跟蹤技術(shù)非冰川區(qū)均RMSE

Tab.4RMSEinice-freeregionofexternalDEM-
assistedoffsettrackingtechniquesusingdifferentincidenceCOSMO-SkyMeddata

(m)

數(shù)據(jù)對(duì)方位向偏移量RMSE距離向偏移量RMSE20160630/201607310.062 0.051

通過對(duì)比TerraSAR-X 2種方法的非冰川區(qū)結(jié)果RMSE可知，當(dāng)SAR影像對(duì)垂直基線較小時(shí)，2種方法的RMSE都比較?。?并且，基于外部DEM輔助的方法RMSE更低，說明在小空間基線時(shí)，DEM輔助的方法也有助于提高結(jié)果的精度。當(dāng)空間基線較大時(shí)，DEM輔助的方法在距離向上改善效果明顯，方位向上也有提升。由公式(1)和(2)可知，主要是因?yàn)榫嚯x向上地形引起的偏移量誤差受到空間基線的影響較大，方位向上受到空間基線的影響較小。依據(jù)表4可知，對(duì)不同入射角的COSMO-SkyMed數(shù)據(jù)應(yīng)用本文方法時(shí)，監(jiān)測(cè)結(jié)果在非冰川區(qū)的RMSE也很小，能獲得精度較高的偏移量結(jié)果。以上分析表明，大空間基線下在地形起伏較大區(qū)域使用DEM輔助的方法能夠有效降低地形引起的偏移量誤差，減小空間基線長(zhǎng)度的限制，提高偏移量跟蹤結(jié)果的精度，擴(kuò)大SAR數(shù)據(jù)的使用范圍。

4 結(jié)論

1)偏移量跟蹤技術(shù)不受冰川區(qū)域影像復(fù)相干性降低的影響，是一種監(jiān)測(cè)山地冰川表面運(yùn)動(dòng)很有效的工具。

2)傳統(tǒng)的偏移量跟蹤技術(shù)忽略了地形起伏對(duì)結(jié)果的影響，在地形變化劇烈以及影像空間基線較大的情況下距離向結(jié)果的精度明顯降低。

3)DEM輔助的偏移量跟蹤技術(shù)在地形陡峭以及SAR數(shù)據(jù)對(duì)空間基線較大時(shí)能夠獲得可靠的結(jié)果，減小了地形以及空間基線長(zhǎng)度帶來的限制，有利于更好開展山地冰川表面運(yùn)動(dòng)監(jiān)測(cè)。

4)不同波束模式的SAR影像應(yīng)用DEM輔助的偏移量跟蹤技術(shù)時(shí)仍能獲得較好的冰川表面運(yùn)動(dòng)監(jiān)測(cè)結(jié)果，進(jìn)一步提高了SAR影像的利用率。

志謝： 衷心感謝GAMMA遙感公司Urs Wegmuller博士耐心的答疑和幫助。