劉 冰,張永紅,吳宏安,康永輝,姜德才
(1. 山東科技大學,山東 青島 266590; 2. 中國測繪科學研究院,北京 100830)
隨著我國經濟發展和城市化進程的加快,高速公路的建設及運營與國民經濟發展水平愈顯密切。截至2015年底,全國高速公路總里程超過12萬千米,基本實現了《國家高速公路網規劃》中的由7條首都放射線、9條南北縱向線和18條東西橫向線組成的“7918”國家高速公路網計劃[1-2]。
快速發展的經濟使得地下水和其他地下資源如石油、天然氣和煤等過度開采,從而導致了地面沉降日益嚴重。以北京市為例,其地面沉降的主要原因是地下水資源的過度開發。北京市是以地下水為主要供水資源的城市,有60%以上水源來源于地下水[3]。地面沉降不僅會造成建筑物的毀壞、地下管道破裂、海水倒灌,還會嚴重威脅城市交通網絡的安全運行,高速公路上的不均勻沉降更會破壞路基結構,長期累積易引發橋頭跳車等事故,嚴重危害人民的生命和財產安全。因此,高速公路網的監測對確保高速公路的安全運營與形變預測具有重要意義。
精密水準觀測和GPS觀測是高速公路沉降監測的常規方法,但其存在覆蓋空白區,作業周期長,且需要耗費大量人力、物力和財力的缺點,難以快速、客觀地獲取大區域地面沉降。1998年Garbriel首次提出了利用DInSAR技術探測地表微小形變,因其具有不受時間空間影響、廣覆蓋、成本低、分辨率高等優點,在地表形變監測中得到了廣泛應用[4],2006年Power[5]等利用DInSAR技術監測了美國3條高速公路的形變情況,2010年文獻[6]將InSAR與GPS技術相融合,獲取了河南省禹登高速采空區的地表形變。但DInSAR技術易受時間空間失相干及大氣相位干擾的影響[7],為克服其缺點,時間序列InSAR技術應運而生。時間序列InSAR技術不同于常規的DInSAR技術,它是選取在較長時間尺度內仍保持高穩定性、高信噪比特性和較強反射性的點目標,通過這些穩定點目標的相位信息,獲取大區域長時間范圍內的地表形變特征[4,7-8]。2012年文獻[9]利用時序SAR技術獲取了江蘇一采空區上方高速公路的形變;2014年文獻[10]通過PS-InSAR技術獲取了京滬高速的地面沉降情況。但常用的時間序列InSAR技術仍有不足,如PS-InSAR[11-13]技術所需影像數目較多(一般要超過25景影像),且受相干性影響較大;SBAS InSAR[14-15]技術未考慮大氣影響;CT-InSAR[16-17]技術雖結合了二者的特點,但非線性參數估計需同時考慮低、高分辨率的影像,計算復雜且容易產生新的誤差[4,18-19]。鑒于以上原因,張永紅于2014年提出了多主影像相干目標小基線干涉技術(multiple images coherent targets small baselines InSAR,MCTSB-InSAR)[20-21]。本次試驗利用MCTSB-InSAR方法開展北京市高速公路網的沉降監測。
本次試驗研究區位于北京市的東南部平原地區,中心地理坐標約為39.9°N,116.6°E,覆蓋順義區、朝陽區、昌平區、海淀區、東城區、西城區、豐臺區、大興區和通州區,如圖1所示,12條高速公路線和一條快速公路線途經此區域。

圖1 研究區范圍及水準點分布位置
本次試驗數據選取德國的TerraSAR-X衛星條帶模式(StripMap,SM)SAR影像,中心頻率為9.65 GHz,數據格式為單視復數據(SLC)。研究區范圍共使用兩景影像覆蓋,面積約為56 km×32 km,其中第1景(左側)影像包含23期數據,時間跨度為2011年6月—2015年6月,雷達入射角約為33.1°,第2景(右側)影像包含25期影像,時間跨度為2011年6月—2015年7月,雷達入射角約為35.3°。兩景影像的時間、空間基線參數見表1。為降低數據運算量,避免數據冗余,本次試驗對整幅影像作了多視處理,多視處理后兩景影像方位向和距離向像元大小約為5.6 m、2.7 m。為消除地形對干涉相位的影響,需輔助DEM數據模擬地形相位,因此,下載了研究區約90 m分辨率的SRTM DEM數據。另外,從北京相關部門獲取了研究區2011—2015年17個水準點的觀測數據用于精度評價,水準點分布位置如圖1所示。

表1 Terra_SAR影像時空基線
MCTSB-InSAR技術不使用單一影像為主影像,而是選擇滿足一定時間、空間基線的SAR影像組合成干涉對,從而降低了對數據量的要求。在高相干目標提取方面,MCTSB-InSAR技術不單一地使用幅度離差閾值法或相干系數閾值法,而是依據數據情況,綜合使用幅度離差、相干系數和平均幅度的“三閾值法”提取高相干目標。在形變參數估計方面,MCTSB-InSAR技術一般將形變分為線性形變和非線性形變,若形變具有強烈的非線性,則可將形變分解為高階多項式形變加上殘余形變[22]。
MCTSB-InSAR技術的主要流程為:基線分析、影像配準、地理編碼、小基線干涉對組合、生成小基線干涉圖、高相干目標提取、形變參數估計等。本文重點研究3個方面的內容。
對于在一定時間內獲取的同一區域的N幅SAR影像,在一定時間基線距和空間基線距的限制條件下,可得到M對干涉對組合,其中N與M的關系為
N/2≤M≤NN-1 /2
(1)
對本次試驗,依據實際情況對第1景和第2景影像均設置時間基線距小于300 d、空間基線距小于300 m的限制條件,檢查所獲取的干涉對的質量,剔除干涉圖相干性較差的影像,最終第1景和第2景獲取的干涉對個數分別為63個、60個。
高相干目標包括點目標和分布式目標,高相干點目標是指高穩定性、高信噪比特性和強反射性的散射體,高相干分布式目標是指在一段時間內仍保持較高相干性的散射體集合。通常利用幅度離差法提取穩定點目標,但這種方法一般需要足夠多的影像(大于30幅),而且該方法以高信噪比像元為基礎,對一些低信噪比像元(如陰影、水體)易產生誤選,因此需借助平均幅度閾值的設置將后向散射性較弱的偽相干點剔除,從而保證選點的可靠性。高相干分布式目標可通過平均相干系數閾值提取,通過設置合適的閾值,將平均相干系數大于該閾值的像元認為是高相干目標。通常情況下,利用幅度離差閾值法提取的點目標也具有高相干性,因此,合理使用幅度離差、平均相干系數、平均幅度“三閾值”法可大大提高選點精確性,減少噪聲點的存在。本次試驗覆蓋范圍為北京主城區,具有較多的建筑物,因此相干性較高,依據實際情況,第1景和第2景最終獲取點目標個數分別為224萬和195萬。
在高相干點目標上構建Delaunay三角網連接相鄰像元,為避免Delaunay三角網方法帶來的巨大計算量,本文使用了一種多層級、不同步長的子網最近鄰點目標快速連接方法,在網絡連接階段采用局部Delaunay三角網連接,通過對影像分塊大小和重疊度的控制(子塊越小,重疊度越高,則子塊數目越多),自由調整連接密度;在剔除低質量連接邊之后,通過多層級子網擴展連接實現網絡的快速連通。局部Delaunay三角網較原方法在保證點目標連接穩健的同時,也降低了網絡邊的冗余,相較于復雜網絡連接方法,該方法的數據處理時間僅為其1/3[23-24]。對Delaunay三角網中相鄰兩像元(設為m,n)之間可建立如下相位差模型
Δφm,nTi=Δδm,nTi+Δβm,n+Δαm,n+Δnm,n
(2)
(3)
式中,Δφm,n為相鄰兩像元相位差;Ti為第i幅干涉圖的時間基線;Δνm,n與Δhm,n分別為線性形變速率和高程誤差增量;Δβm,n為非線性形變相位差;Δαm,n為大氣相位差;Δnm,n為噪聲相位差。考慮到大氣相位在空間上具有較強的相關性,可通過設置大氣相關距離消除大氣相位的影響,一般情況下,大氣相關距離設置為1~3 km時,可認為大氣相位差為零。
為求得線性形變速率和高程誤差增量,使用整體相位相干系數,建立最優化方程
(4)
對Delaunay三角網上所有邊完成最大化求解后,將γ≥0.7的邊作為可靠的連接關系,利用某一參考點的線性形變速度和高程誤差,通過區域增長的方法對Δνm,n與Δhm,n積分,得到各高相干點上線性形變速速率和高程誤差的絕對量。
從差分相位中減去線性形變相位與高程誤差相位后,可得到包括非線性形變相位φnon_linear、大氣相位φatm和噪聲相位φnoise的殘余相位φres,對φres進行時空頻譜特征分析,見表2。

表2 不同相位成分時空頻譜特征
由表2可知,在時間頻域上作低通濾波,可提取非線性形變相位,進而通過最小二乘法及干涉組合關系求得非線性形變量,同時對低通濾波后的殘余相位可通過適當的大氣相關距離窗口(3 km)作空間低通濾波處理,提取大氣相位。
本次試驗基于Kaiser窗函數對φres進行卷積濾波[17],提取非線性形變相位
φnon_linear=φres·h(Kaiser)
(5)
式中,h(Kaiser)為時間低通濾波核。將提取到的非線性形變量與線性形變量相加即可得到完整形變量。
利用上述方法,提取了整個試驗區2011—2015年的地表形變信息,年平均沉降速率如圖2所示。試驗區域為北京市主城區,分布著較多建筑物,因此提取的相干點較為密集。首先分析整個試驗區的地表形變分布,由圖2可以看出,北京市2011—2015年平均沉降速率最大值為-182 mm/a,位于朝陽區黑莊戶鄉(A),地表沉降主要分布于朝陽區東部、通州區西北部、海淀區北部、昌平區東南部、順義區及大興區部分地區。其中,朝陽區東部和通州區西北部連成的大面積沉降漏斗是北京市地表沉降最為嚴重的地區。但不同的沉降中心各年平均沉降速率有所不同,各沉降中心位置分布如圖2所示,表3為這些沉降中心不同時間段的沉降速率變化。

圖2 北京市2011—2015年平均沉降速率及各沉降中心分布
表32011—2015年北京市主要沉降中心各年份沉降速率對比mm/a

沉降中心年份2011—20122012—20132013—20142014—2015黑莊戶鄉(A)-180-184-188-175金盞鄉(B)-167-166-173-165高辛莊(C)-110-111-111-104天通苑(D)-115-117-118-111八仙莊(E)-98-99-99-91上莊鎮(F)-92-92-97-92姚各莊(G)-85-85-87-85首都國際機場(H)-81-84-85-79北野場村(I)-48-49-50-49
從表3可知,北京市2011—2014年沉降速率雖呈逐年增長的趨勢,但增長較為緩慢,這可能與近年來北京市采取的地下水禁采措施有關,而2014—2015年各沉降中心較2013—2014年明顯減緩,尤其是黑戶莊鄉沉降速率由-188 mm/a減小至-175 mm/a。據調查,2014—2015年沉降減緩的主要原因為南水北調中線工程的正式通水。2014年底南水北調中線一期工程向北京供水,使得地下水開采量減少,有效緩解了北京市的水資源危機,截至2015年12月,北京市地下水位平均回升7.5 m。
北京高速公路是國家高速公路網的重要組成部分,本研究區包含了京藏高速、京承高速、機場北線高速、京平高速、首都機場高速、機場第二高速、通燕高速、京哈高速、京石高速、京開高速、京滬高速、京津高速及京通快速公路。根據高速公路的分布,在SAR幅度影像上跟蹤出這些線路的中心位置,并設置30 m的緩沖區,提取緩沖區范圍內的高相干點及沉降信息,圖3為提取的12條高速公路與一條快速公路的地面沉降信息,表4為各高速公路沿線高相干點密度統計表。

圖3 研究區內各高速公路的平均沉降速率

序號名稱長度/km高相干點個數高相干點密度/(個/km)1京藏高速25.52494982京承高速251013413機場北線高速11.333402964京平高速36.794662585首都機場高速18.61292706機場第二高速11.528192457京通快速公路18.430101648通燕高速16.633322019京哈高速39.924716210京石高速4.03709311京開高速24.916156512京滬高速23.36042613京津高速24.64585186
本次試驗中使用的Terra_SAR雷達衛星軌道狀態為升軌、右視,當衛星掃描方向與路面平行時,高相干點密度較高,因此東西向道路的高相干點密度高于南北向道路,而對于南北向的機場第二高速公路,因道路兩側分布較多金屬隔音墻,與路面形成二面角強散射,故其高相干點密度較大。實地調查發現,高速公路上提取到的高相干點主要是道路兩側及中間位置的護欄與路燈、信息牌、廣告牌、高速公路收費站等,在位于城市區域內的高速公路路段,過街天橋也會形成很好的高相干點目標。
受朝陽區東部和通州區西北部連成的大面積沉降漏斗的影響,機場第二高速、京通快速公路、通燕高速及京哈高速與京津高速的西段沉降較為嚴重,其中京哈高速最大沉降速率達到-153 mm/a,為北京市高速公路網沉降最大值;京通快速公路、通燕高速、機場第二高速和京津高速最大分別可達-138、-97、-128、-122;機場北線高速、京平高速、首都機場高速也均有較明顯沉降,沉降速率最大值為-86 mm/a,對位于北京路段的京石高速、京開高速及京藏高速南段則均無明顯沉降現象。
另外,從圖3中也可看出,位于沉降區域內的每條高速公路沿線沉降速率均不同,不均勻沉降較為嚴重。高速公路路基的不均勻沉降會使路基本身產生損壞,進而導致路面發生不均勻形變。路面不均勻沉降超過一定標準后,就會損壞路面結構,減少公路使用壽命,同時不平整的路面會增加行車阻力,降低行車舒適度,影響行車安全。國際上將路面平整度作為評價路面質量及行車舒適的重要指標,我國規定高速公路路面平整度標準差應小于3.5 mm[25]。本文在提取的高速公路高相干點信息的基礎上,以2014—2015年的機場第二高速為例,通過ArcGIS及Matlab軟件獲取了其路面平整度標準差大于3.5 mm的路段,如圖4所示。
據統計,機場第二高速全長11.5 km,其中標準差超過3.5 mm的受損路段為2 km。高速公路路面遭到損壞會嚴重影響行車平穩性,極易引發交通事故,因此路面平整度的獲取對公路行車安全具有重要意義。相比于平整度檢測儀,利用InSAR技術可以大大減少人力物力,更加方便快捷的提取路面平整度標準差,因此InSAR技術在獲取路面平整度信息方面具有良好的應用前景。
為了驗證MCTSB-InSAR技術反演地表沉降信息的可靠性,本文利用從北京市相關部門獲取的17個水準點的觀測數據進行精度評價分析。獲取的水準點數據時間為2011年12月—2015年12月,本文所使用的雷達數據時間跨度為2011年6月—2015年6月,時間上較為吻合。因利用MCTSB-InSAR技術獲取的點目標與水準點的空間位置不可能完全一致,為獲得客觀準確的結果,這里使用鄰近點原則,即以水準點為中心,選距其80 m范圍內的最近一高相干點目標參與精度評價分析。比較水準點與最鄰近點目標年平均沉降速率,結果見表5,二者標準差為4.6 mm/a,驗證了利用MCTSB-InSAR技術獲取地表形變的可靠性。

圖4 2014—2015年機場第二高速及受損路段(右側)平均沉降速率

點名水準測量值InSAR形變速率差值(1)81A-42.9-45.5-2.6(1)83新-3.3-5.4-2.1(3)92A-13.8-13.80(5)150.1-2.8-2.9(6)40-42.6-34.28.4I京通10-122.3-115.37.0I京通2-1.6-5.9-4.3I京西13-10.9-11.7-0.8I京西4新-0.4-3.8-3.4I沙懷3基-0.5-2.0-1.5I沙京13-5.4-1.24.2I順京13-3.0-6.5-3.5II黃清11-0.4-1.8-1.4II黃清9-12.4-6.85.6II蘇舊1-0.6-1.2-0.6II孫通3-111.8-101.310.5II小順2-53.9-56.4-2.5標準差4.6
本文利用MCTSB-InSAR技術獲取了2011—2015年北京市主城區372萬個點的地表形變信息,并分析了各時段的地表沉降時空變化,發現2014—2015年北京市平均沉降速率明顯減緩,這與北京市禁止開采地下水與南水北調中線工程通水密不可分。通過建立空間緩沖區,并依據高速公路上固有的高相干點源提取了北京市高速公路網的沉降信息,發現機場第二高速、京通快速公路、通燕高速及京哈高速與京津高速的西段沉降較為嚴重,機場北線高速、京平高速、首都機場高速也有較明顯沉降,高速公路整體呈現不均勻沉降。以機場第二高速為例,利用InSAR技術獲取了其2014—2015年路面平整度信息,為高速公路養護及判斷行車安全打下基礎。通過對比水準點數據與整個研究區形變量結果,發現MCTSB-InSAR技術獲取的地表沉降信息精度較高。
北京市地表沉降的主因是地下水的過度開采,此外,地表高層建筑的集中建設、重要交通網絡的施工及運營,也使得地表載荷增加,從而造成局部差異性沉降。不均勻地表沉降會對高速公路網等基礎設施造成威脅,嚴重危害其運營環境。因此,北京市在日后的規劃管理中應嚴格控制地下水的開采,同時加大對交通網絡的監測力度,提前對可能出現的危險進行規避。
[1] 董雷宏.“全國一張網”下高速公路運營養護管理新思路探討[C]∥中國公路學會養護與管理分會第七屆學術年會論文集. 西安: 中國公路學會養護與管理分會, 2017: 188-191.
[2] 中國地質調查局. 地質災害調查與監測技術方法論文集[M]. 北京: 中國大地出版社, 2005: 2.
[3] 陳蓓蓓, 宮輝力, 李小娟, 等. 北京地下水系統演化與地面沉降過程[J]. 吉林大學學報(地球科學版), 2012, 42(S1): 373-379.
[4] 張紅, 王超, 吳濤, 等. 基于相干目標的DInSAR方法研究[M]. 北京: 科學出版社, 2009: 4-5.
[5] POWER D, YOUDEN J, ENGLISH J, et al. InSAR Evalua-tion of Landslides in Support of Roadway Design and Realignment[C]∥Proceedings of 2006 IEEE International Conference on Geoscience and Remote Sensing Symposium. Danver, CO, USA: Institute of Electrical and Electronic Engineers, 2006: 3848-3851.
[6] 芮勇勤, 陳佳藝, 丁曉利. 基于InSAR與GPS技術的公路采空區變形監測[J]. 東北大學學報(自然科學版), 2010, 31(12): 1773-1776.
[7] ZEBKER H A, VILLASENOR J. Decorrelation in Inter-ferometric Radar Echoes[J]. IEEE Transactions on Geoscience and Remote Sensing, 1992, 30(5): 950-959.
[8] 李英會. 基于時間序列高分辨率SAR影像的地表形變監測技術研究[D]. 阜新: 遼寧工程技術大學, 2012.
[9] 范洪冬, 鄧喀中, 祝傳廣, 等. 基于時序SAR技術的采空區上方高速公路變形監測及預測方法[J]. 煤炭學報, 2012, 37(11): 1841-1846.
[10] 張學東, 葛大慶, 肖斌, 等. 多軌道集成PS-InSAR監測高速公路沿線地面沉降研究——以京滬高速公路(北京-河北)為例[J]. 測繪通報, 2014(10): 67-69, 104.
[11] FERRETTI A, PRATI C, ROCCA F. Permanent Scatterers in SAR Interferometry[J]. IEEE Transactions on Geoscience and Remote Sensing, 2001, 39(1): 8-20.
[12] FERRETTI A, PRATI C, ROCCA F. Nonlinear Subsid-ence Rate Estimation Using Permanent Scatterers in Differential SAR Interferometry[J]. IEEE Transactions on Geoscience and Remote Sensing, 2000, 38(5): 2202-2212.
[13] KIM J S, KIM D J, KIM S W, et al. Monitoring of Urban Land Surface Subsidence Using PSInSAR[J]. Geosciences Journal, 2007, 11(1): 59-73.
[14] BERARDINO P, FORNARO G, LANARI R, et al. 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.
[15] CASU F, MANZO M, PEPE A, et al. SBAS-DInSAR Analysis of Very Extended Areas: First Results on a 60000-km2 Test Site[J]. IEEE Geoscience and Remote Sensing Letters, 2008, 5(3): 438-442.
[16] MORA O, MALLORQUI J, BROQUETAS A. Linear and Nonlinear Terrain Deformation Maps from a Reduced Set of Interferometric SAR Images[J]. IEEE Transactions on Geoscience and Remote Sensing, 2003, 41(10): 2243-2253.
[17] WU T, WANG C, ZHANG H, et al. Deformation Retrieval in Large Areas Based on Multibaseline DInSAR Algorithm: A Case Study in Cangzhou, Northern China[J]. Interna-tional Journal of Remote Sensing, 2008, 29(12): 3633-3655.
[18] 吳宏安, 張永紅, 陳曉勇, 等. 基于小基線DInSAR技術監測太原市2003~2009年地表形變場[J]. 地球物理學報, 2011, 54(3): 673-680.
[19] 葛大慶, 殷躍平, 王艷, 等. 地面沉降-回彈及地下水位波動的InSAR長時序監測——以德州市為例[J]. 國土資源遙感, 2014, 26(1): 103-109.
[20] 張永紅, 吳宏安, 張利民, 等. 基于高分辨率時間序列衛星SAR影像的交通網絡沉降監測[J]. 大地測量與地球動力學進展(第二輯), 2014: 619-631.
[21] 姜德才, 張永紅, 張繼賢, 等. 天津市地鐵線不均勻地表沉降InSAR監測[J]. 遙感信息, 2017(6).
[22] 張永紅, 吳宏安, 孫廣通. 時間序列InSAR技術中的形變模型研究[J].測繪學報,2012,41(6):864-869, 876.
[23] 吳宏安, 張永紅, 康永輝,等.一種面向時間序列InSAR的不連通子網快速連接方法[J].測繪學報,2016,45(10):1192-1199.
[24] 張永紅, 吳宏安, 康永輝. 京津冀地區1992—2014年三階段地面沉降InSAR監測[J].測繪學報,2016,45(9):1050-1058.
[25] 趙巖. 基于高速公路舒適性的不均勻沉降標準研究[J].武漢理工大學學報(交通科學與工程版),2011,35(6):1245-1247.