楊 斌,劉國祥,姚京川,李 政,馮海龍,毛文飛,張 瑞,吳明杰
(1.中國國家鐵路集團有限公司 工程管理中心,北京 100844; 2.西南交通大學 地球科學與環(huán)境工程學院,四川 成都 610031; 3.中國鐵道科學研究院集團有限公司 鐵道建筑研究所,北京 100081; 4.雄安高速鐵路有限公司,河北 雄安 071800)
地表沉降是自然或人類活動引起的一種具有緩變性、不均勻及不可逆性的地質災害[1-2]。長期以來,因我國北方干旱少雨持續(xù)超采地下水,致使多地存在顯著的不均勻地表沉降。隨著時間的累積,已對鐵路基礎設施結構穩(wěn)定性及軌道平整性等造成嚴重影響[1-3]。據(jù)統(tǒng)計,截至2009年,我國有8萬km2的地區(qū)累計沉降量超過200 mm,包括多個知名城市及周邊區(qū)域。因此,及時把握區(qū)域性沉降的空間分布、及時探明其發(fā)展趨勢與演化態(tài)勢,對于保障重大交通基礎設施安全運營和制定相關管控政策,均具有重要的參考價值和意義。
在地表沉降的監(jiān)測領域,基于時域建模分析的合成孔徑雷達時序干涉(TS-InSAR)技術克服了時空失相干、軌道和地形累積誤差、相位解纏誤差以及大氣延遲誤差等負面因素[4-11],能夠準確地提取長時間尺度上的地表形變信息,現(xiàn)已廣泛運用于各種地災調查、凍土監(jiān)測[12]及地表沉降監(jiān)測等。其中,在算法上最有代表性的技術方法包括Ferretti等于2000年前后提出的永久散射體雷達干涉(Permanent Scatterers InSAR, 即PSInSAR)[13-16],Berardino等[16-17]提出的小基線集(Small Baseline Subsets, SBAS)差分雷達干涉方法,以及Hooper等[18]充分考慮地物的穩(wěn)定性并利用三維解纏法有效消除干涉相位模糊度提出的基于后驗估計的PSI方法[19-21]。而側重于干涉點目標分析的技術方法,則主要包括Werner等[22-23]提出的IPTA算法、Zhang等[24]提出的時域相干點目標分析(TCP-InSAR)等,這些方法均可統(tǒng)稱為PSI方法[25]。需要注意的是,PSI方法在進行PS點探測時,振幅離差指數(shù)(Amplitude Dispersion Index, ADI)小于0.25的像元方可認為對應于PS點目標,且該方法至少需要30景以上 SAR 影像進行PS探測以及形變建模和提取。
為了探究京津冀地區(qū)的不均勻沉降對高速鐵路的影響,本文選取2015年11月至2018年3月間覆蓋研究區(qū)域的51景C波段SAR影像,對影像進行精密配準后[1],利用干涉點目標時序分析方法提取該區(qū)域的沉降信息,并結合該區(qū)域人類活動信息和水文地質資料等對高鐵沿線的沉降漏斗進行演化態(tài)勢分析,為該區(qū)域鐵路安全運營及維護等提供數(shù)據(jù)支撐。
PSI技術以某一影像為主影像,將其他多幅影像配準到該主影像上,并經(jīng)重采樣處理,運用統(tǒng)計分析的方法獲取時序的影像強度和相位信息,并識別出不受時空基線和大氣延遲影響的點目標,其主要提取后向散射特性穩(wěn)定的相干目標,對這些相干目標進行干涉處理,在處理過程中將差分干涉相位φdiff分解為[1]
φdiff=φl+φtopo+φres
(1)
其中,φres=φbaseline+φatm+φnoise+φnl
式中:φl為線性形變差分干涉相位;φtopo為外部DEM引入的地形誤差相位;φres為殘余相位;φbaseline為基線誤差造成的殘余相位;φatm為大氣引起的相位;φnoise為噪聲引起的相位;φnl為非線性引起的相位。
干涉點目標時序分析方法以差分干涉相位φdiff為因變量,高程改正Δhe和線性形變速率Δv為自變量建立相位函數(shù)模型,殘余相位為常量建立的相位回歸分析模型[1]為
φdiff=aΔvT+bΔhe+φres
(2)
式中:T為時間基線;B為空間垂直基線;λ,R和θ分別為雷達波長、雷達至地面的斜距和入射角。
利用式(2)在時空域進行回歸分析,首先解算出Δhe和Δv,然后更新差分干涉相位,去除殘余相位,接著進行迭代回歸分析,最終得到與點目標觀測值最佳擬合的形變速率和沉降值。
干涉點目標時序分析方法主要步驟包括:
(1) 獲取N(N≥30)幅SAR影像數(shù)據(jù)[9],裁剪出研究區(qū)域,選取其中1幅為主影像,將其他影像與主影像進行精密配準;
(2) 生成N-1幅干涉圖,經(jīng)過去除地形相位和參考橢球相位等處理,生成差分干涉圖;
(3) 選取振幅離差指數(shù)(ADI)閾值,選出高振幅的像素作為PS候選點,但考慮到大氣相位(Atmospheric Phase Screen, APS)對形變速率和高程誤差估算模型的影響,選點后,需要估算形變速率、高程誤差以及APS[10];
(4) 根據(jù)Ferretti等提出的二維周期圖(Two Dimensions Periodogram, TDP)[11]解算方法,獲取每個點的線性形變速率和高程誤差;
(5) 對殘余相位時空濾波,分離出大氣相位、非線性形變和噪聲相位(忽略不計),最終得到研究區(qū)域內PS點的形變速率和高程信息。
1.2.1 PS點探測
振幅離差指數(shù)PS點選取是Ferretti等[10]提出,其主要是利用PS點的強反射特性,挑選出高振幅值的像素作為PS的候選點,而劉國祥等[26]通過對PS散射特性的穩(wěn)定性進行分析,提出利用振幅離差指數(shù)閾值從PS候選點中進一步精選出PS。本文利用后者提出的方法進行PS點選取,其選取條件為
(3)
式(3)中第1個公式確定相干性,第2個公式確定的是穩(wěn)定性。
1.2.2 時序建模及線性參數(shù)求解
永久散射體時序差分雷達干涉(PSI)是一種單一主影像的時序InSAR技術,利用1.2.1節(jié)中的方法進行PS點選取,借助外部DEM生成PS點在時間維度上的差分干涉相位。差分相位包括形變相位、大氣延遲相位、高程誤差及隨機噪聲,具體可表達為
φdiff=φref+φtopo+φdef+φatm+φnoise
(4)
式中:φref為參考橢球相位;φtopo為地形相位;φdef為真實形變相位。
為了克服大氣延遲等影響,PSI技術采用鄰域差分建模方法。相鄰2個PS點在第i幅干涉圖中的相位差Δφdiff,i可表示為
(5)
式中:Bi和Ti分別為第i幅干涉對的空間垂直基線和時間基線;(xl,yl)和(xp,yp)分別為第i幅干涉圖中點l和點p的位置坐標; Δε(xl,yl;xp,yp), Δφres,i(xl,yl;xp,yp;Ti)和Δv(xl,yl;xp,yp)分別為第i幅干涉圖中點l和點p的高程誤差、殘留相位增量和LOS向的形變速率增量。
利用Ferretti等[11]提出的二位周期圖解算方法,在|Δφres,i|<π的條件下,可基于多個纏繞的基線增量,通過使目標函數(shù)最大化,即式(6),并搜索估計出二維解空間最優(yōu)的Δε和Δv。
(6)
式中:γ為基線模型相干系數(shù);Δωi為觀測值與擬合值之差。
1.2.3 非線性形變及大氣延遲相位
在差分干涉相位中,除了形變速率和高程誤差外,還有非線性形變相位,大氣延遲相位和噪聲相位,即殘留相位φres,i為
φres,i=φnl,i+φatm,i+φnoise,i
(7)
式中:φnl,i為第i幅干涉圖的非線性形變相位;φatm,i為第i幅干涉圖的大氣延遲相位;φnoise,i為第i幅干涉圖的噪聲。
由于各幅干涉圖相位的殘留相位在時空域的特性不同,因此可對以上幾個相位進行分離。利用空間低通濾波消除干涉失相關和其他隨機相位分量,然后在時間域濾波,便可以通過時序變化趨勢識別對應的非線性形變相位。得到的真實形變相位φdef為
φdef=φl+φnl
(8)
研究區(qū)域為京津冀地區(qū),主要包含雄安縣、霸州市、廊坊市、北京以及部分縣鎮(zhèn)等,面積約為1.5萬(km)2。該研究區(qū)域是京津冀一體化的核心區(qū)域,人口密度較大,經(jīng)濟開發(fā)區(qū)眾多,交通基礎設施較多且路網(wǎng)復雜。區(qū)域內大部分地區(qū)地下水開采處于超采狀態(tài);地質條件復雜,存在地裂縫、崩滑流和地面塌陷、地表沉降以及地下水污染等環(huán)境地質問題[27];地基主要為軟土、松軟土層,壓縮性較高[1]。研究區(qū)域如圖1所示。
圖1 研究區(qū)域
選取研究區(qū)域2015年11月至2018年3月51景C波段的SAR影像,具體干涉對及其時空基線見表1,分辨率(距離向×方位向)約為2.3 m×13.9 m,寬幅約為250 km,2016年12月14日的影像被選為主影像,從影像與主影像的空間基線最長為257.104 6 m,最短為3.291 6 m,而時間基線最長為468 d,最短為12 d。此外,為了去除參考橢球和地形相位的影響,選取30 m空間分辨率的AW3D30數(shù)字地表模型(DSM)數(shù)據(jù)作為地形數(shù)據(jù),該數(shù)據(jù)官方公布的高程精度為5 m標準差。
根據(jù)上面提出的技術方法,得到研究區(qū)域的沉降速率場,如圖2所示。不難發(fā)現(xiàn),研究區(qū)域局部沉降較為明顯,沉降漏斗主要分布在雄縣、霸州市以東的勝芳鎮(zhèn)和左各莊鎮(zhèn)、固安縣、廊坊市以及北京通州區(qū),全局沉降速率分布在20~206 mm·a-1之間。其中沉降最為嚴重的是勝芳鎮(zhèn),最高沉降速率達206 mm·a-1,漏斗中心累積沉降量為248 mm;左各莊鎮(zhèn)、通州地區(qū)、雄縣和廊坊市的沉降漏斗中心沉降速率分別為159,152,185和110 mm·a-1,累計沉降量分為186,181,200和124 mm;而沉降最為緩慢的是固安縣,但沉降速率也達79 mm·a-1,漏斗中心累積沉降量為116 mm。李廣宇等[1]利用2015年6月至2016年8月間的C波段SAR數(shù)據(jù),監(jiān)測到左各莊鎮(zhèn)、通州區(qū)黑莊戶、勝芳鎮(zhèn)和廊坊市的沉降速率分別為163,159,197和108 mm·a-1,累積沉降量分別為191,187,234和120 mm,與本文監(jiān)測結果基本相當。與李廣宇等[1]文中采用的數(shù)據(jù)集相比,本文觀測的時間跨度為2016年11月至2018年3月,新增了2016年9月至2018年3月的數(shù)據(jù),一方面,驗證了前人在研究區(qū)域監(jiān)測結果的可靠性;另一方面,在一定程度上表明2016年9月至2018年3月間,研究區(qū)域的沉降變化趨勢基本穩(wěn)定。由圖2可知,高鐵線1沿線經(jīng)過雄縣、霸州及固安縣等區(qū)域,高鐵線2沿線河北段經(jīng)過了廊坊市,高鐵線3沿線北京段經(jīng)過通州區(qū),而以上分析表明,這幾個區(qū)域均存在年沉降速率超過100 mm·a-1的大范圍沉降漏斗,其中漏斗中心累積沉降量最大的是雄縣,達200 mm。這些不均沉降可能會對高速鐵路的安全運營及相關基礎設施的穩(wěn)定造成潛在的危害,應引起關注。
為了進一步分析高速鐵路沿線區(qū)域的不均勻沉降演變過程,以高鐵線1為例,對鐵路左右兩邊10 km范圍內的沉降進行分析。圖3給出了高鐵線1沿線區(qū)域的年沉降速率。從圖3可以看出,線路沿線的沉降漏斗主要集中在雄縣至黃村段,沉降明顯,范圍分布較大,其中,雄縣至霸州段不均勻沉降最為顯著,漏斗中心相對的最大沉降速率達185 mm·a-1。為了分析鐵路沿線的不均勻沉降變化,在鐵路沿線選取9個PS點(圖3中點a—k) ,對其進行時序分析。
圖3 高鐵線1沿線10 km沉降速率緩沖圖
圖4給出了9個PS點累積形變量的變化情況。由圖4可以看出:從2015年11月8日至2018年3月27日,雄縣至新機場段總體呈下沉趨勢,沉降速率相對穩(wěn)定;點a—f的累積沉降量均超過130 mm,其中點b的累積沉降量達200 mm,點g和h處的沉降速率緩慢,累積沉降量分別為38和68 mm;點a—k的累積形變量曲線中有4個明顯的加速段 (如圖4中粉紅色框所示),分別是2016年4月12日至2016年5月30日(圖4中粉紅框1所示) 、2017年5月31日至2017年6月12日(圖4中粉紅框2所示) 、2017年10月22日至2017年11月03日(圖4中粉紅框3所示)以及2018年1月2日至2018年2月7日(圖4中粉紅框4所示)。而這9個PS點中,除點h外,其余各點基本上都選在鐵路附近,因此,這些不均勻沉降可能會對鐵路安全運營及其基礎配套設施的穩(wěn)定造成潛在的危害,必要時需提前預防。
圖4 沉降區(qū)域內所選PS點的累積沉降曲線
為了探究高鐵線1沿線的沉降速率分布,圖5給出了鐵路沿線沉降速率。由圖5可以看出:整體上,鐵路沿線呈下降趨勢;里程在0~15 km之間(李營至黃村附近)的沉降速率基本在6 mm·a-1左右波動;里程在15~55 km之間(黃村至固安縣附近)的沉降速率出現(xiàn)3次陡崖式下降,位置分別在約DK23,DK30以及DK44處;里程55~85 km的沉降速率穩(wěn)定在12~38 mm·a-1之間,里程85~95 km之間,沉降速率急劇下降,在DK98前后,出現(xiàn)陡崖式下降,沉降速率約達88 mm·a-1,而相應的地理位置是霸州至雄縣段。通過距離量算發(fā)現(xiàn),該結果與圖3和圖4的結果及地理位置基本對應。
圖5 高鐵線1沿線剖面沉降速率
為了更好地解釋高鐵線1沿線的沉降演變過程,還獲取了該區(qū)域的水文地質資料。通過水文地質資料發(fā)現(xiàn),固安及其附近區(qū)域的淺層水位埋深大致分布在15~25 m之間,霸州地區(qū)大致在28~33 m之間,而雄縣地區(qū)大致在13~35 m之間。整個研究區(qū)域地下水埋深變化較大,水位季節(jié)性變化幅度在3~5 m之間,局部區(qū)域達7~9 m,而這些水位變化與不均勻沉降息息相關。表2統(tǒng)計了研究區(qū)域部分的集中供水井。由表2可以看出,供水井的供水面積普遍較廣,深水井多,年供水量均在10萬m3以上。據(jù)調查,李營至北京新機場路段兩側200 m范圍內分布12 個集中供水井,59 個灌溉井,新機場至雄安路段兩側200 m范圍內分布36個集中供水井,125 個灌溉井,由于淺部地層更松散,工程地質性質較深部地層差,多屬中等壓縮性地層,地層結構是黏性土層與砂層交互成層,為雙面排水的地層結構,在相同的水位下降條件下,同樣厚度的地層,淺部地層引起的局部不均勻沉降比深部地層更嚴重一些。通過以上資料分析,一定程度上解釋了雄縣至黃村段不勻沉降嚴重的原因。
圖6給出了高鐵線1—3沿線及周邊的幾個典型的沉降漏斗,并結合沉降漏斗周邊的人口、交通、工業(yè)以及自然資源等對沉降漏斗演化進行分析。固安縣位于北京市、廊坊市、天津市以及涿州市等中間,京九鐵路縱穿全縣,還有106國道、大廣高速等橫貫東西,北面有大興機場,工業(yè)園區(qū)眾多,流動人口較大,常住人口達52萬左右,且該地區(qū)長期的天然氣、地熱以及石油等自然資源開采,是該地區(qū)沉降的主要原因之一。廊坊市的沉降主要集中在廣陽和安次2個區(qū)。已有研究表明,2個區(qū)人口密度大,產業(yè)和水資源豐富[1,28],地熱和地下水的過度開采,導致該區(qū)域持續(xù)的地面沉降;勝芳和左各莊2鎮(zhèn)工業(yè)園區(qū)眾多,區(qū)域人口密度大,人類相關活動頻繁,結合唐嘉等[26]研究,初步認為該因素為霸州市發(fā)生沉降的主要原因之一;雄縣城鎮(zhèn)人口達17.9萬,城鎮(zhèn)化率達47.82%,降水量為600 mm左右,天然氣、石油、地熱及地熱水等儲量豐富,該縣原油年產70萬t,天然氣1 800萬m3,人口的集中,自然資源的大量開采,是該區(qū)域沉降的主要原因之一。前已述及,9各PS點的累積沉降量均存在4個明顯的加速段,而從這4個加速段的季節(jié)分布看,均在夏季或冬季,而夏季用水量巨大,冬季對地熱和天然氣等資源消耗較大,這在一定程度上呼應了以上分析。
表2 研究區(qū)域部分集中供水井統(tǒng)計表
圖6 典型的沉降漏斗
通過對研究區(qū)域全局的沉降速率進行分析,發(fā)現(xiàn)嚴重的沉降漏斗主要分布在雄縣、霸州市、固安縣、廊坊市及北京通州區(qū),最大沉降速率達206 mm·a-1。將本文的監(jiān)測結果與前人相同區(qū)域的監(jiān)測結果進行對比,發(fā)現(xiàn)沉降速率基本相當,表明研究區(qū)域的沉降態(tài)勢趨于穩(wěn)定。分析研究區(qū)域內3條重要的高速鐵路沿線的沉降信息,發(fā)現(xiàn)它們經(jīng)過了雄縣、霸州市、廊坊市及北京通州區(qū)等這幾個沉降嚴重的區(qū)域,漏斗中心的相對沉降速率均在100 mm·a-1以上,且累積沉降量最高達200 mm。利用高鐵線1的沉降信息,發(fā)現(xiàn)雄縣至霸州段不均勻沉降最為顯著,漏斗中心相對的最大沉降速率達185 mm·a-1;分析雄縣至新機場段的9個PS點的累積沉降速率發(fā)現(xiàn),該段整體呈下降勢態(tài),且沉降速率相對穩(wěn)定,但雄縣至霸州段,累積沉降均超過130 mm,且各點在相同的時段出現(xiàn)4明顯的加速沉降段,而4個加速段均在夏季或者冬季;從鐵路沿線的剖面沉降速率來看,其結果基本與沿線的沉降信息及地理位置相對應。為了更好地解釋沉降漏斗的演化態(tài)勢,利用水文地質和相關的人類活動資料以及前人的相關解釋對研究區(qū)域進行分析,結果表明各地區(qū)的不均勻沉降主要與區(qū)域本身的水文地質條件、自然資源的過度開采、區(qū)域人口高度集中、工農業(yè)活動頻繁及交通情況相關。
以上研究和分析表明,研究區(qū)域內的不均勻沉降可能會對高速鐵路的安全運營及相關附屬建筑的穩(wěn)定造成潛在的危害,需要密切關注。時序InSAR技術克服了傳統(tǒng)InSAR技術的不足,能夠精確,并全覆蓋地探測出高速鐵路沿線的沉降漏斗,這是傳統(tǒng)監(jiān)測技術無法相比的,但由于該技術不能進行定位,需借助外部DEM進行地理編碼,完成雷達坐標到地理坐標的轉換,其精度與DEM的分辨率密切相關;另外,星載SAR衛(wèi)星是側視成像,這勢必會產生畸變誤差,影響監(jiān)測精度。因此,聯(lián)合多源數(shù)據(jù)對高速鐵路沿線進行監(jiān)測,是一種必然趨勢。