利用L波段星載重復軌道干涉SAR提取DEM及大氣效應分析

薛東劍，鄭 潔，李成繞，李婉秋，蘇璐璐

(1.成都理工大學 地球科學學院，四川 成都 610059；2.中國林業科學研究院 資源信息研究所,北京 100091；3.國土資源部地學空間信息技術重點實驗室,四川 成都 610059)

利用L波段星載重復軌道干涉SAR提取DEM及大氣效應分析

薛東劍1,2,3，鄭 潔1,3，李成繞1,3，李婉秋1,3，蘇璐璐1,3

(1.成都理工大學 地球科學學院，四川 成都 610059；2.中國林業科學研究院 資源信息研究所,北京 100091；3.國土資源部地學空間信息技術重點實驗室,四川 成都 610059)

以四川省茂縣、安縣地區的ALOS PALSAR L波段雷達影像為數據源，在對地形復雜區的配準、解纏等算法研究的基礎上，結合SAR特有的成像幾何結構，對兩幅SAR圖像進行快速自動配準，配準誤差小于0.2個像元；在去除平地引起的干涉相位變化后，運用MCF對纏繞相位進行解纏，在此基礎上提取研究區的數字高程模型(DEM)，分析發現精度受多種因素影響，其中波長較長的L波段數據比波長較短的時間相關性好,此外為消除大氣效應，采用改進的相位累積法去除大氣的影響，在一定程度上提高DEM精度。其誤差范圍在±10 m左右，對快速提取地形信息具有一定的借鑒意義。

PALSAR；InSAR；配準；相位解纏；大氣效應

以往運用遙感技術提取地形信息，多采用光學圖像，通過立體相對進行提取；圖像畸變及天氣條件嚴重限制了其應用。到20世紀60年代末，合成孔徑雷達干涉測量技術(InSAR)應運而生，其能夠全天候、全天時工作，精確度高，具有光學圖像無法比擬的優勢，有著廣闊的應用前景。Graham(1974)等[1]利用機載SAR數據獲取滿足1∶25 萬地形圖要求的高程數據，開創InSAR技術在對地觀測中獲取三維信息的先河。隨著星載合成孔徑雷達的不斷升空，產生大量豐富的SAR數據，雷達干涉被廣泛地應用于提取地形信息[2-5]。研究表明, 波長較長的L波段SAR 數據比波長較短的X、C波段數據, 在干涉測量時, 時間相關性好, 更適合于地形起伏大、植被覆蓋密、形變梯度大的區域[3]，尤其有利于西南地區高植被覆蓋區地形提取。并且ALOS PALSAR FBS數據的干涉臨界基線距較大，能夠使用更多的影像參與干涉處理。由于SAR特有的成像幾何結構，其處理中的很多關鍵技術還一直處于不斷研究與發展中，尤其針對不同的研究區，其影響提取精度的處理是目前研究的焦點與熱點(如解纏、大氣效應等)。采用茂縣地區的L波段ALOS PALSAR數據，在對配準、解纏算法等研究的基礎上，提取了數字高程模型(DEM)，并重點分析了大氣效應對提取精度的影響。采用該技術可快速提取地形信息，這對于地形變化監測及由地形引起的災害監測預報、環境治理等問題具有重要的實用價值。

1 研究區概況及數據預處理

研究區位于四川省東北部的茂縣、安縣地區，地處龍門山中南段與四川盆地相結合區域，地理坐標：東經103°35′47″～104°22′40″，北緯 31°23′48″～32°00′13″，是“5.12”地震中受災較為嚴重地帶，對其震后地形的提取具有一定的應用價值。研究區地形主要為平壩、丘陵(臺地)、中低山三種類型，既以大光包斜沖斷層和北川沖斷層為界，西北部屬四川西部地槽區的前龍門山褶斷帶，系龍門山脈，地勢較高，最高峰位于九頂山主峰高4 989 m，為區內最高點，最低點位于安縣秀水鎮，海拔為548 m。

研究中采用2010-01-24與2010-03-11 1.1級L波段(波長λ：0.236 m)ALOS PALSAR數據，每景數據覆蓋范圍約為2 565.2 km2。主輔單視復數(SLC)影像距離向像元大小4.68 m ，方位向3.18 m，9 344列，18 432行(見表1)，其中以2010-01-24數據為基準(升軌)，進行干涉處理(見圖1，為研究區ALOS PALSAR HH數據，為顯示方便，在此進行了多視處理及地理編碼)；參考DEM數據選擇美國NASA和NIMA聯合測量的SRTM數據。

圖1 研究區ALOS PALSAR數據示意

成像時間中心經緯度/(°)分辨率/m緯度經度距離向方位向入射角/(°)近距點斜距/km極化方式時間基線/d空間基線/m2010012431．6943648103．97570492010031131．6944510103．9831671468318387402387474846717HH46479

2 InSAR關鍵處理技術

利用星載SAR重復軌道方式進行干涉測量，主要是通過重復兩次觀測，獲得同一區域的單視復數影像對，根據兩幅天線和觀測目標之間的幾何關系(見圖2)，結合觀測平臺的軌道參數，在對相位差解算的基礎上，提取地面高程。

圖2 雷達干涉測量原理示意圖

由圖中幾何關系可推算高程：

h=H-ρcosθ.

(1)

(ρ+Δρ)2=ρ2+B2-2ρBcos(90-θ+α).

(2)

(3)

(4)

圖2中，A1，A2為天線在兩個成像時刻的位置；H是A1的高度；B是基線，α為基線與水平方向的夾角；θ為天線在A1處的視角，λ為波長；ρ，ρ+Δρ為A1,A2到觀測點Z的斜距；h為目標點高程。利用重復軌道獲取的相位差僅是相位差值中的小數部分，需聯合相位解纏求出φt。由于SAR特有的成像方式，使其在數據獲取上具有一定要求，且處理上也有一定的難度，尤其是對兩幅圖像的配準及相位解纏的好壞直接關系到提取地形的精度。

2.1 SAR圖像自動配準

對不同天線獲得數據進行精確配準，估測兩幅SLC影像之間距離和方位向的偏差，是利用SAR進行干涉測量提取DEM的關鍵步驟。為了準確計算對應的地面目標的干涉相位或形變相位，需要實現雷達圖像之間亞像元級精度的配準[12]。在光學圖像上，自動配準的基本思想是相同的地物具有相似的光譜特征，通過對兩個圖像做相對移動，找出其相似性量度值最大或差別最小的位置作為圖像配準的位置。SAR圖像的配準也是類似的原理，目前國內外大部分配準方法基本上是在確定控制的基礎上，以控制點坐標擬合多項式：

(5)

其中，xr,yr為參考影像中控制點坐標；Xt,Yt為待配準數據相應控制點坐標。如何精準確定控制點，是配準處理的前提，也是關鍵，文中針對ALOS PALSAR成像特征，綜合分析相位差平均強度梯度法、相干系數法、絕對差值法、最大頻譜法的基礎上，采用相干系數法進行了控制點選取，此方法充分利用了振幅與相位信息，在一定的搜索范圍內，逐一計算每一點的相關系數，以最大值為配準點的位置。為了避免產生模糊問題，取得精確的偏差估計，在多視估計后，進行一次單視處理，并采用改進的最小二乘法生成配準多項式(表2)，從而使配準誤差小于0.2個像元，提高干涉相關性。

表2 干涉相對配準擬合參數

2.2 干涉圖生成及相位解纏

對干涉像對進行配準后，SLC1與SLC2影像的數據可表示為

μ1=|μ1|ejφ1，μ2=|μ2|ejφ2.

(6)

通過共軛相乘，可以計算出同名點上的相位差：

|μ2||μ2|e-j(φ2-φ1).

(7)

其中，μ*表示復數的共軛，μint的相位是每一同名點上相應于兩個觀測點的相位差(φ)，所得相位的主值數據就是被纏繞相位數據，以干涉圖形式表示。運用InSAR技術提取DEM的精度是否成功，關鍵是干涉像對的相位相關性的高低。經處理，在本次實驗中，相干系數最大值是0.989，最小值是0.001 4，平均值是0.691 6， 0.4以上相關性占78%，該試驗區干涉像對的整體相干性較好(圖3)。此時生成的干涉圖中是由多個成分組成，需去平處理，把干涉圖中的地球彎曲距離和方位向的相位趨勢去除，去除了平地相位之后，干涉條紋由密集變為稀疏，避免了將平地相位信息作為誤差引入高程和形變信息。從圖4可以看出，去平后干涉圖，在地形比較平坦的地區(圖4右上角地區)干涉條紋已經去除，這與實際地形比較相符。

圖3 相干圖與強度合成

圖4 去平后干涉

相位解纏在干涉處理中具有重要的地位，而方法的選擇會很大程度影響高程的精度。常用的算法有區域增長法，最小費用流量法等，其中，區域增長算法，是根據相位的質量來規劃積分路徑，先選擇高質量的區域進行積分，然后再向低質量的區域，從而得到解纏圖像。這種方法在相位圖質量較高時，解纏結果好，速度快，但無法識別殘差點，會引入相位誤差影響精度。最小費用流法是一種基于網絡規劃的相位解纏方法，其基本思想是最小化解纏后相位的導數與纏繞相位的導數之間的差異，該方法可兼顧精度與處理效率，尤其對范圍較大地區的解纏具有一定的優勢。在此以最小費用流量(MCF)算法設置的閥值作為標準(在此相干閾值設為0.3，圖5)，將所有相干性滿足條件的相位添加到一個集合，在集合中建立三角網，構造對偶圖，用最小費用流法連接各殘差點，求出最小費用流集合，再對相位圖進行積分，根據相位質量從高到低進行解纏。

圖5 相位解纏結果

3 地形信息提取及大氣效應分析

經過相位解纏圖像(圖5)的相位主值已經轉化成了真實的相位，采用式3計算出每個像元的高程，在此基礎上將獲取的高程信息由雷達坐標轉化為地理坐標(圖6)。為了分析提取出的DEM精度，采用SRTM DEM數據做參考，進行比較分析。SRTM數據在平原區精度優于1∶250 000 DME,中誤差為4.921。

圖6 InSAR提取的DEM與SRTM DEM疊合分析

運用InSAR技術提取的DEM與SRTM DEM數據比較發現，整體趨勢比較吻合(圖6)，平坦地區誤差范圍較小，而坡度較大的地區，誤差相對較大，存在部分異常值，但差值大的像元比較少。①研究區大部位于山區，有植被覆蓋的山區受時間失相關的影響，對提取的精度造成一定的影響，而城市等平坦地區受時間失相關的影響較小, 其相關性較好；②雷達干涉測量中，除了時間去相干外，大氣效應也是重要的影響因素，不同時間獲取的圖像受大氣影響產生不同的相位延遲，干涉相對中兩個SAR圖像回波信號的相位可表示為

(8)

(9)

其中，ρ1，ρ2為斜距，Δρ1，Δρ2為雷達電磁波穿過大氣層(電離層與對流層)，受電子濃度、氣壓、溫度、水汽等影響造成的大氣延遲，λ為波長。Δφ為兩次觀測大氣引起的相位差：

(10)

大氣造成的延遲主要發生在電離層與對流層。以往研究表明，大氣影響主要是由于水汽的變化引起的[17-18]，表征大氣水汽的物理量主要為可降水汽含量PWV(Precipitable Water Vapor)，可以轉換為天頂濕度延遲ZWD(Zenith Wet Delay)：

(11)

其中，K為折射常數，ρw為液態水密度，TM為對流層的加權平均溫度，Rv為大氣常數。針對重復軌道干涉測量中高程與相位之間的參量關系，可推導出：

(12)

σh為高程誤差，σδ為天頂濕度延遲誤差，θ為入射角，B⊥為垂直基線。

綜合目前研究現狀[15-20]，InSAR測量中常用的減少大氣影響的方法，除了數據本身因素外，主要有相位累積法，PS技術，基于GPS數據的改正方法以及利用外部水汽產品(如MERIS、FY、MODIS等)校正法等，其中PS技術、GPS校正及外部水汽校正法，雖然普適性較高，但主要適用于地表形變監測中[21]。結合ALOS PALSAR數據特點，且考慮到外部數據獲取的同步性及相位解纏的難易等因素，文中借鑒相位累積法，采用2010年ALOS PALSAR數據生成三對獨立干涉圖進行相位梯度累加分析，相位梯度可表示為

(13)

(14)

為減少噪聲的影響，對干涉圖的相位梯度減去平均值和垂直基線的乘積,得到相位梯度變化

(15)

則大氣引起的相位延遲異常▽φatm為

(16)

▽φvari為相位梯度異常,γi為空間相關系數。通過該相位梯度累積平均一定程度上降低了大氣的影響。大氣改正后的高程差異值明顯低于未經過大氣改成的高程差異值，且高程差異分布相對比較集中，高程差異的標準差由15.65 m提高到了6.02 m，尤其是平坦地區誤差范圍在±10 m以內的置信度由68.78%提高到了93.21%，DEM精度得到提高。

4 結束語

以ALOS PALSAR數據為例，在對InSAR處理流程及算法詳細分析的基礎上，研究了地形復雜區域SAR圖像配準、解纏等核心算法；以最小費用流(MCF)算法為基礎，將相干性滿足閾值條件的相位添加到一個集合，構造對偶圖，用最小費用流法連接各殘差點，對相位進行積分，根據相位質量從高到低進行解纏，并提取了DEM。研究表明利用雷達干涉測量提取地形信息具有自動化程度高、效率較高等特點，但失相關是限制INSAR技術應用的一個較為嚴重的問題，尤其植被覆蓋度大、地形復雜地區，時間去相關明顯；除此外大氣效應也是重要的影響因素，在地形提取中，運用相位梯度累加在一定程度上可以消除大氣的影響，提高DEM的精度； ALOS數據由于波長較長，對地表穿透深度大，干涉性相對較好，有利于西南地區高植被覆蓋區地形提取。

[1] GRAHAM L C.Synthetic Interferometer Radar for Topographic Mapping [A].Proceeding of the IEEE,1974,62(6)[C]. [ s.l.]:[ s.n.],1974.763-770.

[2] GOLDSTEIN R M, ZERBER H A,WERNER C L.Satellite Radar Interferometry :Twodimensional PhaseUnwrapping[J]. Radio Science,1988,23(4):713-720.

[3] RAUCOULES D, COLESANTI C, CARNEC C. Use of SAR interfer ometr y fo r detecting and assessing g ro und subsidence[ J] . C.R. Geo science, 2007( 339): 289-302.

[4] ZEBKER H A, ROSEN P A, Hensley S. Atmospheric Effects in Interferometric Synthetic Aperture Radar Surface Deformation and Topographic Maps [J]. Journal of Geophysical Research,1997,103( B4): 7547-7563.

[5] HSIEH Chia-Sheng,SHIH Tian-Yuan,CHING Jyr,et al. Using differential SAR interferometry to map land subsidence: a case study in the Pingtung Plain of SW Taiwan[J], Nat Hazards (2011) 58:1311-1332.

[6] KOCH A,HEIPKE C.Quality assessment of digital surface models derived from the Shuttle Radar Topography Mission(SRTM)[A].the IEEE 2001 International Geoscience and Remote Sensing Symposium[C].Sydney, Australia:2001.

[7] DEO R,MANICKAM S,RAO Y S,et al.Evaluation of interferometric SAR DEMS generated using TanDEM-X data.2013 IEEE International Geoscience and Remote Sensing Symposium(IGARSS)[C].Melbourne,VIC:IEEE,2013:2079-2082.

[8] 陳爾學,李增元,車學儉.SAR 干涉測量數字高程模型提取與高程誤差校正[J].高技術通訊,2000(7):57-63.

[9] 劉國祥,丁曉利,李志林,等.InSAR建立DEM的試驗研究[J].測繪學報,2001,30(4):337-343.

[10] 吳宏安，張永紅.ALOS/PALSAR InSAR 數學模型研究[J].遙感信息,2011(2):106-109.

[11] 沈強,喬學軍,金銀龍,等.ALOS PALSAR雷達影像InSAR數據處理中的基線和地形誤差分析[J].大地測量與地球動力學,2012,32(2):1-6.

[12] 劉廣,郭華東,范景輝.基于外部DEM的InSAR圖像配準方法研究[J].遙感技術與應用,2008,23(1):72-76.

[13] LEE J S.Digital Image Smoothing and the Sigma Filter.ComputerVision[J]. Graphies and Image Proeessing,1983, 24:255-269.

[14] WANG Qi-jie. Generalized functional model of maximum and minimum detectable deformation gradient for PALSAR interferometry[J].Trans. Nonferrous Met. Soc. China 24(2014) 824-832.

[15] HANSSEN R F．Radar Interferometry：Data Interpretation and Error Analysis[M].Dordreeht：Kluwer Academic Press,2001.

[16] 徐佳,關澤群,何秀鳳.重復軌道雷達干涉測量中的大氣影響及其研究展望[J].國土資源遙感,2007,72(2):1-7.

[17] ZEBKER H A,ROSEN P A,HENSLEY S. Atmospheric effects in interferometric synthetic aperture radar surface deformation and topographic maps[J]. Journal of Geophysical Research,1997,102(4):7547-7563.

[18] 孫廣通,張永紅,吳宏安.合成孔徑雷達干涉測量大氣改正研究綜述[J]. 遙感信息,2011(4):111-116.

[19] SONG Xiaogang,LI Deren,SHAN Xinjia,et al. Correction of atmospheric effect in ASAR interferogram using GPS and MODIS data[J]. Chinese Journal of Geophysics, 2009,52(6):1457-1464.

[20] PUYSSEGUR B,MICHEL R,AVOUAC J P,et al.Tropospheric Phase Delay in Interferometric Synthetic Aperture Radar Estimated from Meterological Model and Multispectral Imagery[J].Joural of Geophysical Research,2007,112:120-135.

[21] 崔喜愛,曾琪明,童慶禧.重軌星載InSAR測量中的大氣校正方法綜述[J].遙感技術與應用,2014,29(1):9-17.

[22] GOLDSTEIN R. Atmospheric limitations to repeat-track radar interferometry,G eophys R. es.L ett., 1995,22(18):2517-2520.

[23] SANDWELL D T, PRICE E J. Phase Gradient Approach to Stacking Interferograms [J].Journal of Geophysical Research, 1998, 103 ( B12 ): 30183-30204.

[24] FERRETTI A, PRATI C, ROCCA F.Multibaseline InSAR DEM Reconstruction:the Wavelet Appro ach[J] .IEEE Transactions on Geoscience and Remote Sensing, 1999, 37 (2):705-715.

[25] Ramon Brcic, Alessandro Parizzi, Michael Eineder,et al. Estimation and compensation of onospheric delay for SAR interferometry[C].IEEE International Geoscience & Remote Sensing Symposium,2010,38(1):2908-2911.

UsingrepeatorbitL-bandinterferometricSARtoextractDEMandanalyzingtheatmosphericeffects

XUE Dongjian1,2,3, ZHENG Jie1,3,LI Chengrao1,3, LI Wanqiu1,3, SU Lulu1,3

(1. College of Earth Sciences of CDUT, Chengdu 610059,China；2. Chinese Academy of Forestry, Beijing 100091, China;3.Key Laboratory of Geoscience Spatial Information Technology of Ministry of Land and Resources,Chengdu 610059,China)

This paper, taking ALOS PALSAR L-band radar images as the data source in the counties of Maoxian and Anxian,Sichuan Province，adopts a coherent coefficient method aiming at two SAR images automatic registration and MCF algorithm at phase unwrapping based on the research of registration, unwrapping algorithms, extracted study area’s Digital Elevation Model(DEM). The key technology of interferometric synthetic aperture radar is focused on the basis of the analysis of the current research status. At the time of interferometry, the SAR data with longer wavelength L-band has better temporal coherence than the shorter wavelength.Then the accuracy and characteristics of DEM data are analyzed, whose most error range is about ± 10 m, and the influence of atmospheric effects in interferometriv aperture radar topographic map is discussed, which can provide the basic data for the research work of subsequent surface deformation.

PALSAR;InSAR; image registration;phase unwrapping; atmospheric effect

2016-11-23

四川省教育廳重點資助項目(16ZA0100)；國土資源部地學空間信息技術重點實驗室開放基金(KLGSIT2013-06)

薛東劍(1977-)，男，博士.

著錄：薛東劍，鄭潔，李成繞，等.利用L波段星載重復軌道干涉SAR提取DEM及大氣效應分析[J].測繪工程，2018，27(1)：5-9，14.

10.19349/j.cnki.issn1006-7949.2018.01.002

TP75

A

1006-7949(2018)01-0005-05

李銘娜]