多頻率InSAR提取沼澤濕地DEM精度對比分析

付波霖，耿仁方，李 穎，高二濤，范冬林

(1.桂林理工大學測繪地理信息學院，廣西 桂林 541004；2.南京信息工程大學遙感與測繪工程學院，江蘇 南京 210044；3.中國科學院東北地理與農(nóng)業(yè)生態(tài)研究所，吉林 長春 130102)

濕地是地球上最重要和最有價值的生態(tài)系統(tǒng)之一，被譽為“地球之腎”，具有蓄洪防旱、凈化水質(zhì)、調(diào)節(jié)氣候、維護生物多樣性等重要的生態(tài)功能與服務價值[1-3]。但是，近些年濕地受人類干擾不斷加劇，在全球范圍內(nèi)呈現(xiàn)消失或退化的趨勢，據(jù)不完全統(tǒng)計，世界上大約57%的濕地已經(jīng)轉(zhuǎn)化或消亡[4]。由于濕地的高程數(shù)據(jù)可以為估算濕地水位、生態(tài)需水量、監(jiān)測濕地地下水位和地表水位變化提供基礎(chǔ)數(shù)據(jù)，因此，及時地更新研究區(qū)DEM對于濕地保護和合理開發(fā)起著至關(guān)重要的作用[5-7]。

濕地多處偏遠且難以到達的區(qū)域，無法使用傳統(tǒng)的DEM測量手段，且更新DEM難度大，費時費力，觀測點的空間密度和時間采樣率較難保障；利用光學遙感立體像對生產(chǎn)濕地DEM，容易受到氣候條件的影響[8-9]；LiDAR系統(tǒng)數(shù)據(jù)源單一，觀測面積受限[10]。與之相比，InSAR技術(shù)不受復雜氣候條件的影響，可以快速、高精度、大區(qū)域地進行對地觀測，現(xiàn)已被廣泛應用于地表變形監(jiān)測、災害監(jiān)測、區(qū)域DEM提取和濕地水位監(jiān)測等領(lǐng)域[11-13]。如文獻[14]利用1993—1996年L-band JERS-1 SAR數(shù)據(jù)干涉測量和監(jiān)測了Florida南部沼澤濕地水位變化，結(jié)果表明雷達干涉測量技術(shù)可以精確估算沼澤濕地水位變化，精度可以達到5～10 cm；文獻[15]基于C-band ENVISAT ASAR和L-band ALOS-1 PALSAR單視復數(shù)數(shù)據(jù)干涉測量和監(jiān)測了黃河三角洲天然濕地水位變化，達到了厘米量級的監(jiān)測精度；文獻[16]利用基于歐空局ERS-1/2的TanDEM數(shù)據(jù)成功提取了香港地區(qū)DEM，統(tǒng)計方差為35.2 m；文獻[17]以高分辨的TerraSAR-X影像為數(shù)據(jù)源提取了伊斯坦布爾城市地區(qū)的DEM，結(jié)果表明高分辨率的TerraSAR-X影像可用于DEM提取，并可應用于要求5～10 m垂直精度的城市地區(qū)。綜上所述，InSAR技術(shù)被國內(nèi)外學者廣泛應用于濕地水位監(jiān)測，并取得了厘米級的精度。但是目前基于InSAR提取DEM的研究主要集中在城市地區(qū)，對于利用InSAR能否提取濕地DEM，精度如何，不同波長SAR數(shù)據(jù)對生成的DEM在精度上是否存在差異的相關(guān)研究，國內(nèi)外學者鮮有涉及。

綜上，本文以黑龍江洪河國家級自然保護區(qū)為研究區(qū)，利用InSAR技術(shù)提取沼澤濕地的DEM。為了探究不同波長InSAR提取沼澤濕地區(qū)域DEM的精度差異，選取3種波長干涉SAR數(shù)據(jù)對，并隨機從1∶10 000地形圖中選取111個點數(shù)據(jù)進行精度驗證，最后對比分析沼澤濕地植被對于不同SAR波長的干涉相干性差異。

1 InSAR基本原理

InSAR提取DEM的基本原理如下：利用具有干涉成像能力的SAR傳感器以重復軌道觀測的方式獲取同一地區(qū)具有一定視角差的兩幅或兩幅以上的單視復數(shù)影像對，對其進行干涉處理得到同一個目標對應的兩個回波信號之間的相位差，并結(jié)合精確的軌道參數(shù)，利用衛(wèi)星軌道與地面目標之間相對的幾何關(guān)系獲取高精度、高分辨的地面高程信息，從而重建地表DEM。如圖1所示。

圖1中，設雷達天線中心S1的高程為H1，地面點P的高程為Z，天線中心S1對目標點P成像時的側(cè)視角為θ1，兩天線中心之間的距離(基線距)為B，基線與水平方向的夾角為α，B∥為基線B在水平視線方向的分量，稱為水平基線，B⊥為基線B在垂直視線方向的分量，稱為垂直基線，R1、R2分別為雷達成像時兩天線中心到目標點P的斜距，ΔR為兩斜距之間的差值，ΔR=R1-R2。

假設在t1和t2時刻，兩幅圖像中來自同一目標的散射特性相同，即隨機相位相同，則此時通過復共軛相乘即可獲得復干涉條紋圖，接收信號的相位只與傳播路徑有關(guān)，可表示為

(1)

則真實的干涉相位Δφ主要取決于SAR后向散射信號的路徑差ΔR

(2)

在三角形S1S2P中，由余弦定理可知

(3)

對于星載雷達系統(tǒng)通常B?R1,再由式(2)可得

(4)

則地面點P的高程Z為

(5)

在不考慮兩側(cè)觀測時目標沿雷達視線向移動引起的相位變化、大氣波動帶來的相位延遲及其他相位噪聲的影響下，在已知星載SAR數(shù)據(jù)的基線參數(shù)和干涉參數(shù)等信息的情況下，可由式(5)確定地面點的高程值。

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

2.1 研究區(qū)概況

研究區(qū)為黑龍江洪河國家級自然保護區(qū)，位于黑龍江省三江平原腹地，同江市與撫遠縣交界處，面積21 835 hm2(如圖2所示)。研究區(qū)為寒溫帶大陸性季風氣候，四季分明，有6個月的冰凍期，年平均氣溫為1.9℃，年降水量為585 mm。洪河自然保護區(qū)是三江平原乃至整個東北原始淡水沼澤濕地全貌的縮影和典型代表，也是我國現(xiàn)有的36塊國際濕地之一。研究區(qū)內(nèi)主要有島狀林植被型、灌叢植被型、草甸植被型和沼澤濕地植被型，且在空間上呈環(huán)狀分布。

2.2 遙感數(shù)據(jù)源

本文分別選取3種不同波長的干涉SAR數(shù)據(jù)對：2016年5月28日和6月8日的X-band TerraSAR HH單視復數(shù)數(shù)據(jù)、2015年10月4日和10月16日C-band Sentinel-1A IW模式的VV單視復數(shù)數(shù)據(jù)、2007年9月7日和10月23日L-band ALOS-1 PALSAR精細模式的HH單視復數(shù)數(shù)據(jù)，以上SAR數(shù)據(jù)均能完整覆蓋洪河濕地自然保護區(qū)，具體SAR衛(wèi)星成像參數(shù)見表1。SRTM 30 m DEM數(shù)據(jù)產(chǎn)品來自美國USGS數(shù)據(jù)服務平臺(http:∥earthexplorer.usgs.gov/)，精度驗證數(shù)據(jù)為保護區(qū)1∶10 000地形圖，來自中國科學院東北地理與農(nóng)業(yè)生態(tài)研究所。

表1 研究區(qū)不同SAR衛(wèi)星成像參數(shù)

3 InSAR提取沼澤濕地DEM技術(shù)流程

InSAR技術(shù)提取地表DEM的主要數(shù)據(jù)處理流程包括：單視復SAR影像對的基線估算、主輔影像配準、干涉圖生成、平地相位剔除、干涉圖自適應濾波及相干性計算、干涉圖相位解纏、軌道精煉和重去平，以及初始相位轉(zhuǎn)高程生成DEM。具體的工作流程如圖3所示。

3.1 單視復SAR影像對基線估算

基線估算是利用InSAR技術(shù)提取地表DEM的前提，其通過計算垂直基線、臨界基線、多普勒位移量、高度模糊數(shù)和形變模糊數(shù)等參數(shù)評價兩幅干涉SAR影像對的質(zhì)量狀況。當兩幅SAR影像成像模式不一致或垂直基線超過了臨界基線的長度，則該影像對相干性丟失，不能進行相干處理[18]。當垂直基線距在臨界基線內(nèi)，通常情況下選擇長基線距的兩幅SAR影像對進行相干處理。表2描述了3種SAR數(shù)據(jù)對基線估計的詳細參數(shù)。

3.2 SAR主輔影像對配準

SAR主輔影像對配準是InSAR干涉處理的基礎(chǔ)。在星載重復軌道SAR主輔影像之間的相對偏移量未知的情況下，一般采用粗配準、像元級配準、子像元配準等計算過程實現(xiàn)主輔影像對的配準。由于單視復SAR影像中既含有強度信息又含有相位信息，可以任意將復影像的相干系數(shù)、強度影像的相關(guān)系數(shù)，以及相位差的平方和最小作為測度進行主輔影像的配準。通常情況下，兩幅影像配準的誤差應盡量在1/10個像元之內(nèi)[19-20]。

3.3 主輔影像對預濾波或多視處理

由于干涉影像對在方位向和距離向均存在一定的頻譜偏移，從而會在生成的干涉條紋圖中引入相位噪聲，降低影像對的相干性。為了高精度地配準主輔影像和提高干涉圖的質(zhì)量，可在方位向和距離向?qū)χ鬏o影像進行預濾波處理。在數(shù)據(jù)處理過程中，可根據(jù)頻譜偏移量的大小選擇是否進行濾波處理。

3.4 干涉圖生成

對配準和多視處理之后主、輔影像對進行復共軛相乘生成干涉條紋圖，干涉條紋圖的稀疏反映了地表的起伏狀況。復SAR影像對共軛相乘獲得干涉條紋圖，圖中干涉相位是以2π為周期循環(huán)記錄，并由于受平地效應的影響，干涉條紋在地勢平坦且起伏和緩的區(qū)域呈現(xiàn)有規(guī)律的分布，具體如圖4所示。

從圖4中可以看出，X-band TerraSAR生成的干涉圖，干涉條紋非常密集并存在大量相位噪聲，C-band Sentinel-1A和L-band PALSAR干涉條紋相對稀疏，同時L-band PALSAR生成的干涉條紋圖對保護區(qū)的描述明顯好于波長更短的C-band Sentinel-1A和X-band TerraSAR。

表2 3種SAR數(shù)據(jù)對的基線估計參數(shù)

3.5 平地相位提出

平地效應是指在干涉條紋圖中，平坦地區(qū)的干涉條紋會出現(xiàn)隨距離和方位向的變化而呈周期性變化的現(xiàn)象。剔除平地相位之后的干涉相位可以更好表征地形相位與參考面的相位差，減少干涉圖的條紋頻率，降低干涉圖濾波和相位解纏的難度。本文利用SRTM 30 m DEM模擬試驗區(qū)的平地相位來消除平地效應。

3.6 干涉圖自適應濾波及相干性計算

地面散射特性、影像配準誤差、系統(tǒng)熱噪聲及數(shù)據(jù)處理過程中的誤差均為干涉圖中相位噪聲的主要來源，相位噪聲降低了干涉圖的質(zhì)量，增加了相位解纏的復雜性，直接影響最終生成DEM的精度。因此，本文利用自適應濾波器盡量多地剔除干涉圖中的相位噪聲，并采用最大似然估算器計算相干系數(shù)，相干系數(shù)可以指示區(qū)域相干程度的高低，研究區(qū)SAR影像自適應Goldstein濾波[21]的干涉條紋圖如5所示。

從圖5中可以看出，經(jīng)過自適應Goldstein濾波處理和相位去平之后，在相干性較好的區(qū)域，干涉條紋圖變得更為稀疏和平滑，對于洪河自然保護區(qū)的表征更為詳細；深水沼澤植被區(qū)受濕地水體和體散射的影響，去相干性較為明顯，干涉圖不連續(xù)。

3.7 干涉相位解纏

復共軛相乘獲得的干涉圖中的相位值為干涉相位的主值，此時干涉相位只能以2π為周期，循環(huán)記錄相位變化，必須對其進行相位解纏處理，確定各個像元之間的真實干涉相位值，才能利用干涉圖獲取地面高程信息。

洪河保護區(qū)復雜的植被群分布結(jié)構(gòu)決定了干涉圖存在很多不連續(xù)的區(qū)域，3D Delaunay最小費用流(minimum cost flow)方法采用Delaunay三角剖面進行相位解纏，更加適用于研究區(qū)干涉圖的相位解纏。因此，本文采用最小費用流方法進行相位解纏，獲取干涉條紋中的真實相位值，并且為了讓設定的解纏閾值盡可能考慮到區(qū)域的整體相干情況，數(shù)值設置為0.12～0.15。

3.8 軌道精煉和重去平

由解纏相位轉(zhuǎn)換到地表高程，需要有精確的基線參數(shù)，必須進行軌道精煉，進一步精確估算干涉測量所需的幾何參數(shù)，并再次去除平地效應。軌道精煉的常用方法有：精確衛(wèi)星軌道輸入二次糾正和基于地面參考點的數(shù)學模擬。本文通過選取地面控制點進行軌道精煉和重去平，盡量選擇植被群落較為均一且相干性較高的區(qū)域，糾正之后的均方根誤差控制在1個像元之內(nèi)，經(jīng)過軌道精煉和重去平處理的解纏相位如圖6所示。

3.9 真實相位轉(zhuǎn)高程生成DEM

獲取各個像元之間的真實干涉相位值之后，可以借助式(5)提取地表的DEM數(shù)據(jù)并進行相應的插值和地理編碼處理。

4 結(jié)果分析與驗證

4.1 DEM精度驗證與結(jié)果分析

為了便于進行精度驗證，在洪河自然保護區(qū)內(nèi)，利用ArcGIS fishnet工具生成250個1 km×1 km格網(wǎng)，利用Create Random Points工具生成111個隨機點。將隨機點、格網(wǎng)與1∶10 000地形圖疊加，使1個格網(wǎng)中只有1個具有精確高程值的隨機點，將這些點作為InSAR提取DEM的精度驗證數(shù)據(jù)并統(tǒng)計分析精度差異，對比結(jié)果如圖7和表3所示。

從圖7和表3中可以看出，InSAR提取的DEM精度受到濕地植被及水體去相干性的影響，深水沼澤植被區(qū)相干性差，進而導致生成的DEM數(shù)據(jù)精度較差，其中基于L-band PALSAR干涉測量提取沼澤濕地DEM的精度高于C-band Sentinel-1A、X-band TerraSAR和SRTM DEM數(shù)據(jù)產(chǎn)品，且與1∶10 000地形圖數(shù)據(jù)吻合度較好，76.58%的高程值差異在3 m以內(nèi)。

表3 InSAR干涉生成DEM與1∶10 000地形圖高程數(shù)值差異分布

數(shù)據(jù)驗證點的高程數(shù)值差<1m<3m<5m<10m<20m<30mPALSAR2485108111111111Sentinel-1A312254692109TerraSAR-X000175111

4.2 不同波長干涉圖中濕地植被相干性差異

解纏相位圖在經(jīng)過軌道精煉和重去平地效應處理之后，計算得到相干系數(shù)圖(如圖8所示)。相干系數(shù)反映了單視復數(shù)SAR影像對所生成干涉圖質(zhì)量的好壞，也可以定量估計SAR兩次成像過程中目標相位的穩(wěn)定性，相干系數(shù)數(shù)值范圍是0～1，數(shù)值越大表明干涉圖的質(zhì)量越好。圖8對應不同干涉SAR影像對生成的相干系數(shù)圖，圖8(a)中剖面線對應的相干系數(shù)值如圖9所示。

從圖8和9中可以看到，在洪河自然保護區(qū)內(nèi)，島狀林和灌草結(jié)合的植被分布區(qū)相干系數(shù)值較大，其次是淺水沼澤植被區(qū)，最小的是深水沼澤植被區(qū)。造成上述相干性差異的主要原因是濕地水位在成像時間內(nèi)的變化導致了以浮水植物為主的深水沼澤區(qū)與部分淺水沼澤區(qū)的植被位置和高度也發(fā)生了變化，而陸生植被除自然生長之外較為穩(wěn)定。從L-band PALSAR和C-band Sentinel-1A相干系數(shù)曲線圖中可以看出島狀林和灌草植被區(qū)相干系數(shù)值達到0.6～0.8，淺水沼澤植被區(qū)相干系數(shù)值為0.4～0.6，深水沼澤植被區(qū)相干系數(shù)值為0.2～0.3。而在X-band TerraSAR提取的相干系數(shù)曲線圖中，島狀林和灌草結(jié)合的植被區(qū)相干系數(shù)值達到0.3～0.4，淺水沼澤植被區(qū)相干系數(shù)值為0.2～0.3，深水沼澤植被區(qū)相干系數(shù)值為0.1～0.2。同一植被區(qū)相干系數(shù)值的差異表明：L-band PALSAR數(shù)據(jù)對更適合沼澤濕地的干涉測量，X-band TerraSAR數(shù)據(jù)對不適合用于沼澤濕地干涉測量。

5 結(jié) 語

本文利用3種不同波長InSAR數(shù)據(jù)對提取沼澤濕地區(qū)域的DEM，使用Goldstein濾波方法消除由時間和基線產(chǎn)生的相位噪聲，提高了干涉條紋的可見性，降低了相位解纏的復雜性；使用3D Delaunay最小費用流方法進行相位解纏，考慮到區(qū)域的整體相干情況，解纏閾值設置為0.12～0.15；選取控制點進行軌道精煉和重去平時，盡量選擇植被群落較為均一且相干性較高的區(qū)域，糾正之后的均方根誤差控制在1個像元之內(nèi)。DEM精度驗證結(jié)果表明：L-band PALSAR更加適合利用雷達干涉測量技術(shù)提取沼澤濕地的DEM，其與1∶10 000地形圖數(shù)據(jù)吻合度較好，76.58%的高程值差異在3 m以內(nèi)，波長更短的C-band Sentinel-1A和X-band TerraSAR由于在時間去相干和濕地植被生長產(chǎn)生的去相干效應，降低了DEM的數(shù)值精度；通過對比同一植被區(qū)相干系數(shù)值的差異，發(fā)現(xiàn)L-band PALSAR相干系數(shù)比其他兩個波段的數(shù)據(jù)要高，更適合沼澤濕地的干涉測量；通過對比不同濕地植被類型的相干系數(shù)，發(fā)現(xiàn)島狀林和灌草結(jié)合的植被分布區(qū)相干系數(shù)值較大，其次是淺水沼澤植被區(qū)，最小的是深水沼澤植被區(qū)。