基于Stacking-InSAR 和SBAS-InSAR 地面沉降監測方法研究
——以新疆阿勒泰地區為例

孫 鵬

（新疆水利水電勘測設計研究院有限責任公司，新疆烏魯木齊 830000）

0 引言

合成孔徑雷達干涉測量技術（Interferometric Synthetic Aperture Radar，InSAR）因其具有高精度、高分辨率、全天候等優點，是監測地表形變中一種方法。該技術應用在城市地表沉降監測、礦區地表沉降監測、山體滑坡監測等項目中，發揮了極大的作用。InSAR 技術演化出Stacking-InSAR 技術與SBAS-InSAR 技術，并可將二者結合起來作為一種組合方法應用在地表沉降、滑坡監測中［1-2］。

本文利用InSAR 技術對新疆阿勒泰地區布爾津滑坡進行監測與分析。Sentinel-1 系統在周邊區域擁有豐富的InSAR 數據積累、高重訪周期與大面積的影像覆蓋范圍，收集基于2014 年10 月20 日至2021 年2 月9 日采集的Sentinel-1 升軌影像，數據處理以增強型SBAS-InSAR 集成時序分析為主要技術手段，獲取布爾津滑坡區域形變的趨勢、范圍和梯度等信息，并從中篩選出地表形變變化量達到預警臨界值，和發生形變突變的可疑區域，作為地表形變監測預警的重點工作區域。

1 基本原理與技術優勢

1.1 Stacking-InSAR技術

干涉圖疊加（Stacking-InSAR）技術的數據處理基礎是進行經典二軌差分InSAR 技術處理。Stacking-InSAR 技術是在獲取M 幅相位解纏的差分干涉圖后，將解纏相位加權疊加，估計出平均相位變化速率。由于大氣相位分量在空間域具有低頻特性、時間域具有高頻特性，故Stacking 技術最大程度地減少了大氣誤差，提高了形變速率解算精度［3-4］。該方法假設每一幅獨立的解纏相位圖中，在不考慮其他噪聲的分布情況下，區域的形變速率為總形變量與時間的比值，即線性形變速率，大氣延遲誤差隨機且相等。通常情況下，為保證結果的準確性，選擇時間和空間基線都比較短的數據進行組合。

1.2 SBAS-InSAR技術

時序分析采用小基線集干涉測量技術（Small Baseline Subset，SBAS）充分考慮大量的干涉圖信息，能有效地消除和削弱解纏粗差、大氣誤差以及DEM誤差等因素的影響，從而獲取更高精度的地表形變成果［5］?；赟BAS 技術的數據處理主要解算步驟如下：

1）對N 幅SAR 影像數據按一定時空基線條件進行干涉組合處理，形成M 幅干涉圖，利用已有DEM（本次采用SRTM）作為外部高程數據，進行差分處理，生成差分干涉圖，以去除地形及平地效應影響；

2）對差分干涉圖進行自適應濾波處理以去除相位噪聲影響，對濾波后的差分干涉圖進行相位解纏，生成解纏相位圖；

3）根據基線條件和干涉相位信息估算高程誤差及線性形變相位，在原始干涉相位中減去估算的線性形變相位得到殘余相位，此時的殘余相位中主要包括非線性形變及大氣相位；

4）解纏此殘余相位，并補償線性形變相位部分得到完整的形變相位，此時相位中還包含有大氣相位的影響；

5）對去除線性形變的殘余相位進行空域低通和時域高通濾波處理以分離出大氣相位；

6）在形變相位中減去大氣相位影響，得到形變相位值；

7）基于SVD 的形變求解；

8）對形變進行地理編碼，獲取WGS-84 坐標系下的形變成果圖。

1.3 總體技術路線

本文根據監測區域基本情況，計劃出合理可行的技術路線。應用雷達監測數據，采用Stacking-InSAR 技術與SBAS-InSAR 技術相結合方式實現了變形數據提取。根據地質變形情況，及時建立災害預警與防治措施。本文總體技術路線如圖1 所示。

圖1 Stacking-InSAR 和SBAS-InSAR 時序分析提取地表形變技術路線圖

2 研究區域與數據源

2.1 研究區域

布爾津滑坡位于新疆阿勒泰地區布爾津縣境內，滑坡位于水電站庫區的上游區域。該區域山體受雨水的沖刷導致地質情況不穩定，很容易出現山體滑坡現象，造成河流堵塞帶來的發電站隱患問題。研究區域范圍內無居住區，地勢起伏較大；最低點海拔+203.4 m，最高點海拔+1125.9 m，平均海拔+841.9 m。該滑坡的研究范圍如圖2 所示。

圖2 研究區域范圍

為比較本文提及的Stacking-InSAR 和SBASInSAR 技術在滑坡體監測中的有效性，在滑坡區域布設監測點，并采用傳統極坐標法獲取監測點位移。為便于體現各監測點位移方向，建立監測坐標系，X方向為平行河道，Y方向為垂直河道，通過坐標變化量直觀表示各點位移趨勢，監測基準網及監測點布置如圖3 所示。

圖3 監測基準網及監測點布置圖

2.2 數據源

2.2.1 影像數據

為避免SAR 影像幾何畸變的影響，所用SAR數據為升軌Sentinel-1 數據，采用TOPS（Terrain Observation With Progressive Scans）掃描方式，共計153 景，覆蓋時間段為2014 年10 月20 日至2021年2 月9 日。數據基本參數如表1 所示，影像采集范圍如圖4 所示。

表1 升軌SAR 數據的基本參數

圖4 研究范圍數據覆蓋圖

2.2.2 DEM獲取

形變監測過程中，數字高程模型（DEM）是非常重要的外部數據，在生成差分干涉圖時需要DEM 來模擬地形相位，此外在雷達坐標系和地理坐標系相互轉化（地理編碼）的過程中也需要高精度的DEM。

美國太空總署（NASA）和國防部國家測繪局（NIMA）于2000 年聯合發射的“奮進號”航天飛機完成了為期11 天的航天飛機雷達制圖任務（Shuttle Radar Topography Mission，SRTM）。經過兩年的數據處理，最終獲得了平面精度±20 m，高程精度±16 m 的全球數字高程模型，分辨率包含1 弧秒和3 弧秒兩種，分別對應30 m 和90 m。本項目選用30 m分辨率的DEM 數據用于地形相位的消除。

2.2.3 極坐標差分法獲取數據

依據153 景的SAR 影像的成像時間，選取2014 年10 月20 日至2021 年2 月9 日之間的SAR數據，每次監測方法保持一致，采用相同觀測人員、儀器使用極坐標差分法觀測。結合每次觀測成果，分析得出各監測點觀測成果相對首期成果水平位移量和垂直位移量。位移量的X方向為正表示沿河向上游方向，反之沿河向下游方向，Y方向為正表示垂直河道向河道方向，反之垂直河道向山頂方向；H方向為0 表示無變化，為負表示下沉。為了驗證Stacking-InSAR 與SBAS-InSAR 對沉降邊緣與沉降中心的監測能力，選取了圖2 中的6 個特征點，如表2 所示。其中，1-1、1-3、2-1 和2-2 屬于沉降邊緣區域，分布于滑坡區域西北側；1-3 和1-4屬于沉降中心區域，2-3 和2-4 位于滑坡區東南側。

表2 極坐標法獲取的變形點坐標值

3 數據處理與結果分析

對獲取的153 景Sentinel-1 數據進行雷達差分干涉測量處理，按照時間序列構建152 個干涉對，其中2020.08.15-2020.08.29 的垂直基線最長，達到了200.89 m，2016.04.15-2016.04.20 的垂直基線最短，數值為43.18 m。按照圖1 所示的Stacking-InSAR 和SBAS-InSAR 技術處理流程對152 個干涉對處理，獲取到相鄰影像之間的形變信息，通過時序累加，得到各時間節點上的累計沉降情況。

依據干涉圖的相干性最大原則，選取2020 年8月20 日的影像為主影像，其余影像與其配準，分別設定時間基線閾值為100 d，空間基線閾值為220 m。將所有滿足閾值的限定的兩景SAR 影像兩兩組合，形成了152 個時序差分干涉對進行處理。按照本文的技術路線圖，采用Stacking-InSAR 與SBAS-InSAR 技術處理152 個干涉對，得到研究區域形變信息。布爾津滑坡153 景SAR 數據計算獲得的2014 年10 月至2021 年2 月地表形變累計沉降量如圖5 所示。布爾津滑坡2014 年10 月至2021年2 月地表形變速率如圖6 所示。從這兩幅圖可以看出：布爾津滑坡區域這段時期內最大累計形變量達-550 mm，最大形變速率達-76.9 mm/a。圖6中，負值表示地表正在遠離衛星，正值表示地表正在靠近衛星。從圖6 可以看出：新疆阿勒泰布爾津滑坡區域分布有三處明顯變形體，變形體1 平均形變-6.0 mm/a，最大形變-25.4 mm/a；變形體2 平均形變-4.8 mm/a，最大形變-17.1 mm/a；變形體3 平均形變-19.0 mm/a，最大形變-76.9 mm/a。

圖5 地表累計形變圖

圖6 2014-2021 平均形變速率圖

為了更直觀感受研究區域內的時空演變規律，這里展示SBAS-InSAR 的時序結果。2014 年10月、12 月與2021 年1 月、2 月時空演變情況如圖7-圖10 所示。

圖7 2014 年10 月時空演變結果

圖8 2014 年12 月時空演變結果

圖9 2021 年1 月時空演變結果

圖10 2021 年2 月時空演變結果

為進一步研究該滑坡形變的時間演化，在3 個形變體上分別提取3 個特征點，獲得滑坡上這些點2014 年10 月至2021 年2 月的時間序列形變，如圖11-圖16 所示。變形體1 最大的累積形變出現在其P3 點上，這個時間段內累積形變達到了-174mm；變形體2 最大的累積形變出現在其P2 點上，這個時間段內，累積形變達到了-223 mm；變形體3 最大的累積形變出現在其P3 點上，這個時間段內，累積形變達到了-499 mm。3 個變形體在2014 年10 月至2021 年2 月存在明顯的非線性形變趨勢?？梢钥吹?，在2018 年1 月至2021 年2月，3 個變形體均出現快速形變現象，形變量比較大，其中變形體3 從2016 年年初就開始出現加速的情況。

圖11 變形體1 特征點示意圖

圖12 變形體1 地表變形速率及特征點累計變形量

圖13 變形體2 特征點示意圖

圖14 變形體2 地表形變速率及特征點累計變形量

圖15 變形體3 特征點示意圖

圖16 變形體3 地表形變速率及特征點累計變形量

4 精度驗證

4.1 傳統測量儀器結果比對

為了驗證Stacking-InSAR 和SBAS-InSAR 技術在滑坡監測方面的精度，采用傳統測繪方法在滑坡體上布置一定的測點?？紤]到滑坡體不同位置移動情況不同，在沉降邊緣區域選擇1-1 號點位，該點位于滑坡體西北側；沉降中心區域選擇1-3、1-4 號點位；2-4 號點位位于滑坡區東南側。各點位分別比較了水平方向與垂直方向移動量。以傳統測繪儀器獲取的數值為真值，為了減小位置差錯產生的誤差，以監測點為中心，將一定范圍內的SBAS點求取平均值作為SBAS 監測值。表3、表4 統計了在2018 年1 月至2021 年2 月時間段內二者比較結果。

表3 水平方向移動量比較

表4 垂直方向移動量比較

由表3、表4 可知：不管是用InSAR 技術還是傳統測繪儀器技術，變形體位移下降趨勢相同，位移量相近、數據相差較小?？梢哉J為采用Stacking-InSAR 和SBAS-InSAR 技術的方式監測地表形變能達到亞厘米級的精度。產生誤差的主要原因是兩種方法監測點位置存在偏差，不能完全對應。采用InSAR 技術監測地面沉降，精度能夠達到厘米級甚至毫米級別，能夠滿足監測項目需求。

5 結語

本文利用2014 年10 月至2021 年2 月期間覆蓋新疆阿勒泰地區的多景Sentinel-1 數據，應用Stacking-InSAR 與SBAS-InSAR 技術相結合的方式獲取了該地區該時間段的地表累計形變量圖，分析得知該地區最大累計形變量達-550 mm，最大形變速率達-76.9 mm/a，可見最大變形體一直處于持續的形變狀態。提出以下建議：

1）實地獲取下游水電站相關數據，對該滑坡地質情況進行調查，做好防護工作與監測風險的情況的工作。

2）為保證變形體能夠得到實時監測，可結合自動化監測技術，保證變形體全方位、多角度、可視化監測。