基于PS-InSAR與SBAS-InSAR的地面沉降監測分析

陳媛媛 朱晨瑋 鄭加柱 王詩涵

摘要：

為了全面掌握長江岸線區域的地面沉降情況，為維護長江防洪安全以及河勢穩定提供參考依據，以長江南京段沿岸5 km緩沖區域為研究區域，采用PS-InSAR和SBAS-InSAR兩種時序分析技術對2017年4月至2021年4月期間30景Sentinel-1A影像進行處理與分析比較，并分析了造成沉降的原因。結果表明：兩種方法所得的特征點的形變曲線變化趨勢一致，通過對沉降點的年平均地表形變速率進行線性擬合以及對其差值進行正態分布分析，驗證了兩種方法所得結果具有一致性和可靠性；兩種方法所得沉降區域分布情況大致吻合，該區域有5個沉降漏斗，分別位于浦口區的長江沿岸區域、雨花臺區東北部與建鄴區西南部的交界處、建鄴區東北部與鼓樓區西南部交界區域、八卦洲區域以及六合區張營水庫區域；地質構造、城市改造建設及承壓水頭下降均對該區域地面沉降產生一定的影響。

關 鍵 詞：

地面沉降； PS-InSAR； SBAS-InSAR； 時序分析； 長江岸線； 長江南京段

中圖法分類號： P237；P642.26

文獻標志碼： A

DOI：10.16232/j.cnki.1001-4179.2023.08.022

0 引 言

江蘇地處長江下游，是南水北調東線工程的水源地，也是長江黃金水道水上運輸大動脈的重要組成部分。江蘇省防洪保安任務重、河勢控制難度大、水資源保護要求高、岸線利用強度大。確保長江江蘇段防洪安全、維護河勢穩定、促進長江生態向好發展等問題是長江治理的關鍵問題之一［1］，其中，確保沿江岸線及其附近區域地面的穩定對于長江治理與保護至關重要，有必要開展長江沿岸區域的長周期持續地面監測工作。

目前，地面沉降監測的手段和技術主要有精密水準測量和全球定位系統（GPS）技術。這些傳統的監測技術存在工作量大、效率低、易受測區范圍影響等問題，且只能獲取離散監測點的下沉量，無法監測地面下沉趨勢和具體范圍。合成孔徑雷達差分干涉測量技術（InSAR）是一種全新的對地觀測技術，它使用兩幅或多幅合成孔徑雷達影像圖，根據衛星或飛機接收到的回波相位差來生成地表形變圖，以探取地表的微小形變信息［2］。InSAR技術具有全天時全天候探測、探測范圍廣、探測精度高、探測成本低等優點，但常規InSAR技術常面臨時間失相關、空間失相關和大氣延遲的問題。為了解決這些問題，就產生了永久性散射體（Persistent Scatterer，PS）和小基線集（Small Baseline Subsets，SBAS）時序分析等技術［3-5］，并且在城市地面沉降、礦山沉降、地震形變、火山活動等領域得到廣泛的應用。目前，已有不少專家基于不同的InSAR技術對南京市的地面沉降情況進行了研究與分析，如：黃其歡等采用13景ERS2 SAR影像，基于短基線D-InSAR方法對南京河西漫灘區進行監測，結果表明該區域最大累計沉降量達12 cm［6］。高二濤等采用2015～2017年覆蓋南京地區的23幅Sentinel-1A影像，基于PS-InSAR和SBAS-InSAR對南京主城區的地面沉降進行監測，結果發現該區域的最大沉降速率達-30 mm/a［7］。楊振等利用SBAS-InSAR技術對2015～2018年的19景南京河西的Sentinel-1A數據進行處理，發現該區域存在4個沉降漏斗［8］。但是，鮮有針對長江岸線南京段區域的時間序列地面監測研究。

本文以長江岸線南京段附近區域為典型研究區域，獲取了30景歐洲航天局（以下簡稱“歐空局”）（Sentinel-1A）衛星雷達影像，基于PS-InSAR和SBAS-InSAR兩種時序分析技術對南京段長江岸線區域進行地面沉降監測，并通過對特征點在兩種方法結果中的形變曲線進行趨勢分析，以及對8 659個經緯度完全相同的特征點的年平均沉降速率及其差值進行相關性分析和正態分布擬合，驗證兩種方法所得結果的一致性，從而保證研究區域沉降監測結果的可靠性與真實性。此外，還從地質構造、城市改造建設等方面對造成長江南京段地面沉降的原因進行了分析。

1 研究區域及數據簡介

1.1 研究區域概況

南京市位于長江下游中部，江蘇省西南部。市域地理坐標為北緯31°14′～32°37′、東經118°22′～119°14′。南京市跨江而居，北邊是廣闊的江淮平原，東邊靠近長江三角洲。《南京市長江岸線保護辦法》對長江岸線做出了解釋：長江岸線指本市行政區域內長江（含洲島）水陸邊界一定范圍內的帶狀區域。本文將長江沿岸5 km范圍內的區域視為長江岸線區域（見圖1紅色邊界內區域）。

南京的地質在全國大地構造單元上屬揚子古陸的北部邊緣，基底主要是輕變質的片巖和變質的火成巖。對于南京長江沿岸地區而言，其地質條件主要特點是地質松散、承載力差、地下水位埋藏淺等。近年來地鐵建設施工及地下空間開發等一系列工程項目的開展，造成南京市地表以及長江岸線區域的沉降日益明顯。

1.2 實驗數據

Sentinel-1A衛星由歐空局于2014年4月在法屬圭亞那庫魯發射，軌道高度693 km，搭載的是C波段合成孔徑雷達。它具有12 d的重訪周期，并且有著多種極化方式和工作模式。本文選用2017年4月至2021年4月南京市30景Sentinel-1A干涉寬幅（IW）模式VV極化的影像進行實驗，影像空間分辨率為5 m×20 m（方位向×距離向）。由于IW模式的影像覆蓋范圍比較大，本文在實驗前對數據進行了裁剪。此外，實驗還獲取了Sentinel-1A衛星對應成像時間的精密軌道數據，地形參考數據使用的是南京地區12.5 m高分辨率DEM產品數據。

2 數據處理

2.1 PS-InSAR技術處理流程

永久散射體的干涉測量技術（PS-InSAR）于2000年由Ferretti等人首次提出［3］，指通過對研究區域內長時間保持穩定散射特征的物體進行研究，即研究一組影像數據中具有高相干性與穩定反射特性的地物點作為PS-InSAR點以提取相位變化，從而避免時間失相關與空間失相關的問題，并且更好地削弱了大氣效應對形變信息的影響，以得到更高精度的地面形變信息與DEM［9-10］。PS-InSAR技術數據處理流程如圖2所示。

本實驗選擇2019年9月2日的SAR影像作為主影像，其余的SAR影像與之組成29對干涉影像對?？臻g基線主要分布于-150～110 m之間，最長空間基線為-140.39 m，最短空間基線為-4.48 m。對每一對像對進行干涉流數據處理，分別為影像配準、干涉圖生成、去平地效應以及振幅離差指數計算。

采用振幅離散指數方法［11］，顯示單個像元值的相干性，輔助PS-InSAR點的篩選，完成高相干性目標點的初步選擇。生成的干涉圖需要剔除偏移相位，偏移相位根據振幅離差指數篩選多個參考點計算得到，參考點的個數可根據研究區域的面積大小設置，本次實驗設置為每25 km2的區域選擇一個參考點，重疊比例參數設置為30%。使用ENVI提供的線性反演模型，基于選擇出的參考點，對地表形變速率與地形殘差進行估算。

根據大氣影響在空間上具有高相關性、在時間上具有低相干性的特性，將反演得到的地表形變速率與地形殘差在空間上采用高通濾波，在時間上采用低通濾波，完成對大氣延遲相位的估算，消除干涉相位中的大氣延遲相位，克服D-InSAR的大氣相干與時空相干等問題。再次使用ENVI提供的線性反演模型重新對地表形變速率與地形殘差進行二次估算，獲得PS-InSAR點最終的形變速率與殘余高程［12］。

為便于在地理空間中對PS-InSAR點進行幾何信息與形變信息的特征分析，需要將上述過程中獲取的雷達坐標系下的PS點進行地理編碼。本次實驗使用12.5 m高分辨率DEM產品對PS點進行地理編碼，相干系數閾值設置為0.75，經過數據處理得到具有高相干性的PS點，在ArcGIS中進行時序分析與空間分析。

2.2 SBAS-InSAR技術處理流程

小基線集算法（SBAS-InSAR）于2002年由Berardino等人提出［4］，該技術僅需要少量的SAR影像數據組合成若干個短基線的干涉相對，從而得到多個差分干涉集合以增加干涉圖的數量，再利用最小二乘法與奇異值分解法在解決時空失相關問題的同時增加采樣率，以此獲取連續且高精度的地表形變信息［13-14］。SBAS-InSAR的主要技術流程如圖3所示。

該部分實驗選擇2019年11月1日的SAR影像作為超級主影像，將其他的SAR影像與其進行配準。配準完成之后，需要對所有的干涉對進行時間基線與空間基線的篩選，由于本次研究獲取的SAR影像數據集時間跨度較長（48個月），故設置時間基線閾值為550 d。此外，因為Sentinel-1A衛星軌道較為穩定，產生的偏差較小，所以設置空間基線為臨界基線的45%。本研究通過引入12.5 m的DEM數據來消除平地效應，采用最小費用流（MCF）方法進行相位解纏［15］。對解纏后的結果進行軌道精煉與重去平處理，選取高質量的地面控制點，通過多項式法來消除軌道參數誤差與去除平地效應。隨后，選擇線性模型來初步估算形變速率與地形誤差，再通過時間域高通濾波與空間域的低通濾波去除大氣延遲誤差。

將SBAS-InSAR反演結果進行地理編碼，將斜距坐標系轉換為地理坐標系。此外，需要將結果的柵格數據轉換為矢量數據，以便在ArcGIS中進行空間分析。

3 結果分析與比較

3.1 PS-InSAR結果分析

通過PS-InSAR數據處理流程得到研究區域的年平均形變速率圖，以長江為對象構建一個5 km的緩沖區對研究區域的地面沉降監測結果進行裁剪，最終得到長江南京段沿線區域的年平均形變速率圖（見圖4），其中紅色矩形表示主要沉降區域，紫色十字符表示特征點位置。

從圖4可以看出，南京市長江流域的年平均形變速率分布范圍為-25～15 mm/a，采用紅色區域表示地面沉降，藍色區域表示地面抬升。研究區域中的主要沉降區域分布于浦口區的長江沿岸區域（A）、雨花臺區東北部與建鄴區西南部的交界處（B）、建鄴區東北部與鼓樓區西南部交界區域（C）、八卦洲區域（D）以及六合區張營水庫區域（E），沉降嚴重的區域大多分布于長江沿岸，最大沉降速率達到-22 mm/a，最大沉降量達到-33 mm，最大累計沉降量為-101 mm。

根據主要沉降區域的分布位置，發現沉降區域A位于浦口區的江浦街道與頂山街道的連接地區，該區域范圍內有地鐵十號線通過且設有3個站點，其最大沉降速率為-19.3 mm/a。包含地鐵二號線與地鐵十號線的建鄴區地面沉降程度最嚴重，其西南地區與東北地區分別有兩個沉降漏斗B與C，且在空間上具有相連成片的趨勢，其中B區域最大沉降速率為-18.2 mm/a，C區域最大沉降速率為-19.8 mm/a。由于長江沖積作用形成的八卦洲地區的南部地區以及位于八卦洲北側六合區的張營水庫附近分別存在大面積沉降區域D與E，其中D區域最大沉降速率為-13.3 mm/a，E區域最大沉降速率為-12.8 mm/a。

3.2 SBAS-InSAR結果分析

圖5為采用SBAS-InSAR技術處理所得的長江南京段緩沖區的形變速率圖，其中紅色矩形表示與圖4對應的沉降區域，紫色十字符表示對應的特征點位置。由該圖可以看出，南京市長江岸線年平均形變速率范圍為-19～10mm/a，在整個長江岸線，有5個主要沉降區，其位置分別在浦口區長江岸線（A）、建鄴區東北部與鼓樓區西南部交界區域（C）、八卦洲區域（D）以及六合區張營水庫區域（E）。最大沉降速率達到-19 mm/a，最大沉降量達到-44 mm，最大累計沉降量為-88 mm，其中，A區（浦口區長江岸線）沉降范圍最廣，C區是建鄴區東北部與鼓樓區西南部交界處的秦淮區莫愁湖區域，其沉降最嚴重。

根據SBAS-InSAR結果中主要沉降區域的分布位置，可以看出A區域是由江浦街道、頂山街道、泰山街道以及沿江街道共4個街道組成的沉降區，該范圍內共有3條地鐵經過且設有7個站點，A區域最大沉降速率為-19.2mm/a；建鄴區存在兩個沉降漏斗B和C，B區域最大沉降速率為-11.9 mm/a，C區域最大沉降速率為-14.4 mm/a；沉降區域D和E的最大沉降速率分別為-12.4 mm/a和-11.8 mm/a。

3.3 兩種技術手段結果對比分析

由于本次實驗缺乏南京市長江岸線同期水準數據，故通過對比分析兩種處理方法的結果來驗證結果的可靠性。選取沉降速率較大的5個區域（A，B，C，D，E區域，見圖4～5），并在每個區域內選取可以代表該區域沉降特征的特征點（圖4～5中的紫色點位），通過繪制5個特征點在整個時間跨度中的累計沉降量變化圖（見圖6）來對比分析兩種技術方法獲取的形變結果。

從圖6可以直觀地看出：特征點的形變曲線整體呈非線性變化，都具有形變緩慢、沉降速率變化不大的特點。PS-InSAR處理方法所得特征點的沉降曲線較為粗糙，形變曲線抖動較大。而SBAS-InSAR處理所得形變曲線更光滑，趨勢也更為穩定。從形變曲線的趨勢和形變結果來看，特征點處的PS-InSAR的最大沉降量為-60 mm，稍小于SBAS-InSAR的-65 mm。兩種方法獲取的結果在沉降量級上雖然存在微小差異，但總體而言具有較強的一致性，沉降區域大致重合，因此實驗結果具有較高的可靠性。

此外，本文還對PS-InSAR與SBAS-InSAR兩種方法提取的8 659個經緯度完全相同的特征點的年平均沉降速率及其差值進行正態分布擬合和相關性分析，來進一步驗證兩種方法獲取結果的可靠性和一致性。如圖7所示，PS-InSAR和SBAS-InSAR兩種方法獲取的年平均沉降速率差值的平均值和標準差分別為0.42 mm/a和0.9 mm/a，滿足CH/T 6006-2018《時間序列InSRA地表形變監測數據處理規范》中的精度要求［16］。圖8中以PS-InSAR監測結果為縱坐標，以SBAS-InSAR監測結果為橫坐標，8 659個沉降點均勻分布在擬合直線兩側，決定系數R=0.789 2，均方根誤差為0.456 1 mm/a，由此可得兩種處理結果具有一致性。

3.4 成因分析

（1） 地質構造。長江岸線南京段兩側漫灘寬廣，據其地形地貌單元與巖土組合，可將其歸為河漫灘軟土地區，主要由長江沖擊、淤積作用形成的第四系沉積軟弱土層構成，故其具有地質松散、承載力差等特點。實驗結果中研究區域有5個主要沉降區，均分布在長江岸線附近，岸線的沉降與地質構造具有密切聯系。其中A區域（浦口區岸線）沉降范圍最大，最大沉降速率達-16 mm，累計沉降量達-40 mm。B區與C區均為三角洲區域，長期受到河水沖擊作用，使得岸底累積大量沉積物，是造成區域沉降的重要因素之一。

（2） 城市改造建設。近年來，南京市地鐵建設及地下工程施工十分頻繁，使得部分區域的沉降日益明顯，其中C區周邊地下有多條地鐵線路運營，該區域內的最大沉降速率達-19 mm/a，累計沉降量達-60 mm。A區域為南京江北新區，作為南京著力打造的長三角區域金融中心核心功能區及國家級發展新區，近幾年建筑改造施工頻繁，地下開發嚴重，而開發區域都是以河灘為主的軟基區域，這是造成該區域大面積沉降的重要原因之一。

（3） 其他因素。相關研究表明，長江岸線南京段地面沉降與該區域附近承壓水頭下降存在較高的相關性［17］。此外，由于其他人為與自然因素，在某種程度上對區域沉降也有著重要的影響。E區為六合區張營水庫，該區域內的地面沉降與壩體的自身結構有關。

4 結 論

本文利用從歐空局獲取的30景Sentinel-1A衛星雷達影像，采用PS-InSAR和SBAS-InSAR時序分析技術對南京長江岸線區域進行地面沉降監測。實驗結果發現：該區域有5個沉降漏斗，分別位于浦口區的長江沿岸區域、雨花臺區東北部與建鄴區西南部的交界處、建鄴區東北部與鼓樓區西南部交界區域、八卦洲區域以及六合區張營水庫區域。對兩種方法所獲結果的特征點沉降速率進行比較，分析了年平均地表形變速率的線性擬合曲線及差值正態分布曲線，發現兩種方法獲取的地面形變結果具有較高的相關性，沉降趨勢與沉降區域的空間分布具有一致性，側面驗證了這兩種方法的可靠性。本文研究可為政府對長江岸線環境保護制定相關決策提供強有力的數據支持，為確保長江岸線區域的防洪安全、維護河勢穩定提供一定的理論依據。

參考文獻：

［1］ 朱宇馳.近50年長江南京段河道變遷研究［D］.南京：南京大學，2016.

［2］ 朱建軍，李志偉，胡俊.InSAR變形監測方法與研究進展［J］.測繪學報，2017，46（10）：1717-1733.

［3］ FERRETTI A，PRATI C，ROCCA，et a1.Nonlinear subsidence rate estimation using Permanent Scatters in differential SAR Interferometry［J］.IEEE Transactions on Geoscience and Remote Sensing，2000，38（5）：2202-2212.

［4］ BERARDINO P，FORNARO G，LANARI R，et a1.A new algorithm for surface deformation monitoring based on small baseline Differential SAR interferograms［J］.IEEE Transactions on Geoscience and Remote Sensing，2002，40（11）：2375-2383.

［5］ 侯安業，張景發，劉斌，等.PS-InSAR與SBAS-InSAR監測地表沉降的比較研究［J］.大地測量與地球動力學，2012，32（4）：125-128，134.

［6］ 黃其歡，徐佳.短基線DInSAR法長江漫灘區大面積沉降監測研究［J］.山東科技大學學報（自然科學版），2011，30（1）：7-10，15.

［7］ 高二濤，范冬林，付波霖，等.基于PS-InSAR和SBAS技術監測南京市地面沉降［J］.大地測量與地球動力學，2019，39（2）：158-163.

［8］ 楊振，徐佳.利用Sentinel-1A數據監測南京河西地區地面沉降［J］.測繪通報，2020（1）：61-65，75.

［9］ REHAN K.Monitoring subsidence in urban area by PSInSAR：A case study of Abbottabad City，Northern Pakistan［J］.Remote Sensing，2021，13（9）：1651-1651.

［10］ 朱邦彥，唐超，任志忠，等.基于PS-InSAR技術的珠海市地表形變監測與驅動力分析［J］.測繪通報，2022（6）：108-113.

［11］ 劉昱彤.PS-InSAR技術的地表沉降監測研究［D］.唐山：華北理工大學，2020.

［12］ 羅三明，杜凱夫，暢柳，等.基于PS-InSAR方法反演北京地區地表沉降速率［J］.大地測量與地球動力學，2014，34（1）：43-46.

［13］ 婁連惠，劉強，譚玉敏.基于高分辨率SAR的三峽庫區地表形變分析：以巴東城區為例［J］.人民長江，2021，52（5）：114-119.

［14］ 高文錦.基于SBAS-InSAR技術的深圳西部地區地表形變研究［D］.西安：西安科技大學，2021.

［15］ 王志勇，孟慶穎.InSAR相位解纏算法在礦區沉降監測中的對比與分析［J］.遙感信息，2014，29（6）：47-51，57.

［16］ 國家測繪地理信息局.時間序列InSAR地表形變監測數據處理規范：CH/T 6006-2018［S］.北京：測繪出版社，2018.

［17］ 朱邦彥，姚馮宇，孫靜雯，等.利用InSAR與地質數據綜合分析南京河西地面沉降的演化特征和成因［J］.武漢大學學報（信息科學版），2020，45（3）：442-450.

（編輯：高小雲）

Abstract：

In order to fully understand the land subsidence in the riparian area along the Changjiang River，and provide references for guaranteeing the flood control safety and the river regime stability of Changjiang River，a 5km buffer area along the Nanjing section of Changjiang River was taken as the study area，two time series analysis technologies，PS InSAR and SBAS InSAR，were used to process，analyze and compare with 30 Sentinel-1A images from April 2017 to April 2021，and the causes of subsidence were analyzed.The results showed that the deformation curves of feature points obtained by the two methods were basically consistent，and the consistency and reliability of the two results were also verified by linear fitting for annual mean ground deformation rate and the normal distribution analysis for the difference of the ground deformation rate.The distribution of the subsidence areas obtained by the two methods coincided roughly.There were 5 subsidence funnels in this area，which were located along the Changjiang River in Pukou District，the junction of the northeast of Yuhuatai District and the southwest of Jianye District，the junction of the northeast of Jianye District and the southwest of Gulou District，Baguazhou and Zhangying reservoir area of Liuhe District.It was found that geologic structures，urban renewal construction and drop of confined water all had an effect on land subsidence in the region.

Key words：

land subsidence；PS-InSAR；SBAS-InSAR；time series analysis；Changjiang River bank；Nanjing section of the Changjiang River